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Heat Transfer Between Concentric 
Rotating Cylinders 


In this paper the problem of convective heat transfer between concentric rotating 
cylinders is studied. Experimental heat-transfer data are presented for four different 
values of clearance between the cylinders and for several combinations of outer to inner 
cylinder speed. The heat-transfer performance indicates three regimes of flow; the first 
at low peripheral velocities, in which laminar flow and heat transfer by conduction 
prevail, the second at cylinder speeds above a theoretically predictable value, in which 
vortex flow occurs and ts the controlling mechanism, and a third at still higher speeds, 
in which a distorted type of vortex motion may be present. The data for the case of the 
inner cylinder only rotating can be correlated by the equation 


Ny u Nyx u 


'/s 
cond 0.175 N Ta 


for the Taylor number range 90 to 2000. A heat-and-momentum-transfer-analogy 
solution for this case follows the trend of the data, but gives results which are somewhat 
high. The combined effects of rotation of both cylinders may be correlated by the em- 
pirical equation 


— — — 3.5 d/R) 
” 
41.1 + (Nn, — d/Ri) 
” 


for values of the abscissa from about 2 to 50. 


Objectives and Introduction 


a A PREVIOUS PAPER [1]! the heat transfer from 
a single rotating cylinder to still air and to a cross-stream of air 
was investigated. Analysis of the thermal behavior of machinery 
often involves not only heat transfer from rotating elements to a 
surrounding fluid of large expanse, but also heat transfer from con- 
fined rotating parts. Such a problem might arise in connection 
with consideration of the thermal behavior of enclosed shafting, 
multiple concentric drives, turbine rotors, and similar machine 


1 Numbers in brackets designate References at end of paper. 

Contributed by the Heat Transfer Division and presented at the 
Annual Meeting, New York, N. Y., November 30-December 5, 1958, 
of THe AMERICAN Society oF MECHANICAL ENGINEERS. 

Note: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those of 
the Society. Manuscript received at ASME Headquarters, August 7, 
1958. Paper No. 58—A-99. 


Nomenclature 


parts. Heat transfer from a long rotating cylinder to a concentric 
cylindrical shell is geometrically the simplest of the possible con- 
figurations, and perhaps one of the most widely encountered; it 
has thus been chosen as the subject for investigation. 

The objectives of the present work have been to study the heat- 
transfer characteristics of such a system and to determine the in- 
fluences of the ratio of outer to inner cylinder angular velocity, 7, 
and of the relative radial clearance between the cylinders, d//,, 
on the heat-transfer rate. Consideration is restricted to the case of 
no axial flow of fluid in the annulus between the cylinders. 

A number of authors have dealt with the flow between concen- 
tric cylinders with the inner one rotating. Taylor [2], in a very 
early paper, showed that above a critical speed of rotation the 
laminar Couette flow breaks down into a flow consisting of pairs 
of three-dimensional ring-shaped vortexes. Hagerty [3] and 
Kaye and Elgar [4] have also photographed such vortex patterns. 
A sketch of this flow pattern is shown in Fig. I(a). 

Taylor [5] also measured the velocity and temperature profiles 


English Letter Symbols 
= surface heat-transfer 
heated cylinder, ft? 


unit 
area of 


convection 
surface of heated evlinder, Btu/- 
(hr ft? deg F) 


conductance on velocity, ft/sec 
distance from surface of heated 


cylinder to a point the 


specific heat at constant pressure, 
Btu/(ib deg F) 

radial clearance between concen- 
tric cylinders, ft 

friction factor for concentric cylin- 
ders, dimensionless 

proportionality factor in Newton’s 
Second Law, 32.2 (Iby)/(Ibr) 
ft/sec? 
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unit thermal conductivity, Btu/- 
(hr ft? deg F/ft) 

heat-transfer rate, Btu /hr 

radius of heated cylinder, ft 

radial co-ordinate in 
laver, ft 

temperature, deg F 

velocity relative to heated evlinder, 
ft /see 


boundary 


boundary layer, ft 


Greek Letter Symbols 
a@ = thermal diffusivity, k/pc,, ft?/hr 
€y = eddy diffusivity for heat, ft?/hr 
€y = eddy diffusivity for momentum, 
ft?/hr 
(Continued on next page) 
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occurring at very high rotational speeds and found that both the 
circumferential velocity and the moment of momentum ap- 
parently increased in moving toward the outer cylinder. He 
supposed that this was due to the effect of the total-pressure tube 
used in the measurements, and by extrapolating to zero tube di- 
ameter he concluded that the true velocity distribution was one 
of constant moment of momentum, or uk = constant, through 
about 80 per cent of the gap. From a comparison with the 
measured temperature profile it was proposed that vorticity 
rather than momentum is transferred by turbulent diffusion in the 
annular space between the two cylinders. 

Much later work by Pai [6] indicates that this may not be the 
case, but that the phenomenon which Taylor observed was quite 
likely due to the persistence of the vortex flow, even at speeds 
several hundred times the critical speed at which it begins. 
Pai’s measurements with a hot wire indicate that a distorted vor- 
tex pattern, as shown in Fig. 1(b), exists and that the mean- 
velocity profile obtained depends upon where the velocity tra- 
verses are made relative to the positions of the vortexes. If 
mean-velocity profiles for several positions along the length of 
the cylinder are averaged, the averaged curve is found to be one 
of constant mean velocity over a large part of the gap. The 
secondary motion is found to play an essential part in the trans- 
port of moment of momentum, and in regions where secondary 
flow is negligible, the momentum follows the laws for homologous 
turbulence. 

Velocity measurements by Hromas and Thompson [7] agree 
with those of Pai, although their measurements with the inner 
cylinder rotating are very limited, the bulk of their data being 
taken in the laminar flow region with the outer cylinder rotating. 

A great deal of attention has been given to the determination 
of the point at which vortex flow begins. Taylor [2] has derived 
an equation for the speed at which laminar flow becomes un- 
stable for small disturbances when the two cylinders are rotated 
in the same direction with a relative velocity. 
pressed in the notation of the present work as, 


This may be ex- 


Rotating Inner Cy! 


Ww Stationery Outer ANN 


(@) Taylor number slightly above critics] value 


‘Rotati Inner Cylinder \ 
Turoulent Turbulent 
Region Region 
Turbulent 


Regton 


Stationary Outer Cyitnder 


(0) Very high value of Taylor number 


Fig. 1 Probable vortex patterns between concentric rotating cylinders 


Meksyn [9, 10, 11], in a series of papers, has presented a number 
of stability criteria. His expression for the critical speed with 
the cylinders rotating in the same direction is somewhat more 
complex than that given by Taylor, but reduces to the same 
value for limiting cases. For the cylinders rotating in opposite 
directions his results may be expressed as: 


Nn 


Via + — 9) — [1 — 


9 


48.7012 + d/R;) 
1+ 0.00056 (1) 
(1 — — + 0.0571 ( 0.652d/R, + 
1+” 
0.652d/R, 


where Ny, is a dimensional group shown by Kaye and Elgar [4] 
to be a significant parameter in concentric rotating cylinder work. 

Taylor's experiments and those of Lewis [8] confirm this ex- 
pression quite well. Lewis’ experiments indicate that it is valid 
for values of d/R, at least as high as 0.71, and is also valid for 
small negative values of » (cylinders rotating in opposite direc- 
tions). The limiting value of critical Taylor number given by 
this expression, for the outer cylinder stationary (n = 
very small clearance ratio (d/R, — 0), is 41.1. 


0) and 


Nomenclature 


This equation is valid provided that (—7) is not much smaller 
than one half. 

Critical values calculated from these expressions are shown with 
the present experimental results in Figs. 6 through 9, where it 
may be seen that they predict quite accurately the transition 
from laminar flow performance. 

Friction data for concentric cylinders with either the inner 
cylinder or the outer cylinder rotating have also been obtained 
by Taylor [12]. Some of his data for the inner cylinder rotating, 


ratio of outer cylinder rpm to inner 


cylinder rpm; 7 is(+)forcoro- ut = 
tation and ( —) for counterrota- 
tion y+ = 

vy = kinematic viscosity, pu/p, ft®/hr 

p = fluid density, lb/ft® 

= fluid shear stress, lb/ft? 

w, = angular velocity of inner cylinder, Nwu 
rad /sec Ner = 

w, = angular velocity of outer cylinder, Nn = 


rad /sec 
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Nondimensional Groupings 
nondimensional relative velocity 
in boundary layer, u/ V 9.7% ‘p 
nondimensional distance from 
heated cylinder surface, 
V gero/p/v 
= Nusselt number, hd/k 
Prandtl number, ye,/k 
Taylor number, 


(Rw,d/v) 


Taylor number for transition 


N Ta, = 


from laminar to vortex flow, 
when cylinder speed ratio has 
the value 7 


Nore: All fluid properties in the ex- 
perimental correlations have been eval- 
uated at a temperature halfway between 
the surface temperatures of the two cylin- 
ders. 
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a 
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Wegrer's Data 
© - a/R, = 0.0308 


4 - a/R, = 0.0524 
Oo - a/R) = 0.0808 
Taylor's Deta 


a/R, = 0.0278 


= 0.0588 


= C.1018 
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Fig. 2 Concentric cylinder friction data of Taylor and of Wagner, inner cylinder only rotating 


along with a large number of data taken by Wagner [13] for 
cylinders with small clearance ratios, are shown on Fig. 2. These 
data have been replotted in terms of a grouping of variables which 
correlates the data for various clearance ratios with a single 
curve. The laminar flow friction factor used in the ordinate may 
be expressed as: 


_ 


(f Nt) tam 


Fig. 2 indicates, in addition to the transition at a Taylor number 
of about 41, what appears to be a second transition at a Taylor 
number of about 90. 

Becker[15], in a very recent work, gives some experimental heat- 
transfer data for concentric cylinders with the inner cylinder only 
rotating.? The only previously available data for this case are 
those which have resulted from experiments related to the cooling 
of electric motors. Most of these latter data have been for slotted 


* Becker’s data were discovered too late to include for comparison 
in Fig. 10. They cover one clearance ratio, d/R; = 0.233, a Taylor 
number range to about 500, and correspond very well with the data 
of this paper. 
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surfaces or for axial air flow in the annulus, but some of the data 
by Gazley [14] are for smooth cylindrical surfaces and no axial 
air flow. These data are included on Fig. 10 for comparison with 
the results of the present investigation and, along with an analogy 
solution proposed by Gazley, will be discussed in a later section 
of this paper. 


Description of Test Apparatus 


The details of construction of the test system are shown on 
Fig. 3. The inner cylinder has been discussed in a previous 
paper [1] and will be described only briefly here. It consists of 
an electrically heated test section 18'/, in. long and 2.26 in. in 
diameter. Plastic spacers insulate the copper section from the 
steel stub shafts at the ends, and the entire inner cylinder as- 
sembly is supported by double sets of ball bearings to minimize 
bending moment at the joints. Heater power is supplied through 
a Variae control to slip rings on the cylinder, and measured by 
means of a high-current watt-meter. Rotation is provided by 
vee-belt drive from a d-c motor on a variable voltage source. 

The temperature of the heated cylinder is sensed by a station- 
ary thermocouple which projects into a cavity in the solid end of 
the rotating copper test section. This system eliminated the 


CHROMALOK CARTRIDGE HEATING ELEMENT 
MODEL C-816, GRASS CASING, SOLDERED IN PLACE 


JRON- CONSTANTAN THERMO - 
——COUPLE, 30GA. DUPLEX 


~ ROTAT NG 
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+ 
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BALSAWOOD THERMOC OU PLE 
HOLDER, NON- ROTATING 
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A- BELT 


Construction details of assembled cylinders 
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need for slip-rings in the thermocouple circuit and was found to 
indicate adequately the local temperature of the copper cylinder. 
The measured temperature is corrected for variation along the 
length of the cylinder due to end leakage to obtain the average 
temperature of the heated cylinder. 

The outer aluminum cylinder is mounted concentrically on the 
inner cylinder by means of duplexed ball bearings, and the axial 
flow of air through the space between the cylinders is obstructed 
by soft leather rings glued to the plastic spacers at the ends of the 
heated-cylinder test section. 

By means of suitable gears and/or pulleys the outer cylinder 
is driven from the same motor as used to rotate the heated cyl- 
inder, By altering the train of gears and/or pulleys, either cylin- 
der can be rotated separately, or both may be rotated simul- 
taneously at fixed speed ratios either in the same or in opposite 
directions. The speed of rotation of each cylinder is measured by 
means of a Strobotac. 

The temperature of the outer cylinder is measured by means 
of a thermocouple imbedded at the point where the local tempera- 
ture is closely equal to the average temperature of the active 
length of the aluminum cylinder, and thus no correction for tem- 
perature variation is necessary. The thermocouple leads are 
brought out to iron and constantan contacts on the surface of the 
cylinder, and iron and constantan pickups are used to complete 
the circuit to the potentiometer. 

It was found necessary to take the reading of this temperature 
with the outer cylinder stationary, since the use of slip-rings with 
the cylinder rotating gave erroneous readings. This was ap- 
parently due to the large diameter, and consequently high surface 
speed, that was necessary for the slip-rings. This measurement is 
made by quickly stopping the cylinder at the end of a run, com- 
pleting the thermocouple circuit by pivoting the pickups into 
place against the contacts on the cylinder, and reading the po- 
tentiometer; the entire procedure takes about half a minute. 
Because of the large time constant of the rather thick-walled 
aluminum cylinder, this method of measurement does not intro- 
duce significant error into the results. An analysis making rather 
severe assumptions indicates that the temperature difference be- 
tween the cylinder and the surroundings can not drop as much 
as | per cent in the time required for the measurement. 

The various gap widths used in the tests were all accomplished 
using the same cylinders. Between each series of runs the 
aluminum outer cylinder was simply dismounted and its inner 
surface rebored to the larger diameter required for the next 
clearance ratio to be examined. 

A crossflow of air to cool the outer cylinder is provided through 
a duct which discharges vertically downward a few inches above 
the test cylinders, 

From an analysis of experimental uncertainty it may be con- 
cluded that the test system is capable of yielding Nusselt number 
determinations within an uncertainty interval of +4'/, per cent. 
The Taylor number determination is within an uncertainty in- 
terval of about +3'/2 per cent. A possible source of error which 
does not lend itself to analysis arises in connection with the large 
time constant of the cylinder, especially at low rotational speeds 
and large clearance ratio. The system does not respond quickly 
to variations in line voltage and variations in surface conductance 
accompanying speed variation. To minimize this source of error, 
a regulated voltage supply was used in the heater circuit, and the 
d-c drive motor was supplied with power from a portable ac-de 
motor generator set which helped to isolate the system from 
variations in the supply. Great care was taken to see that a 
steady state had been achieved before readings were taken. 
Errors arising from unsteady-state conditions would be random 
in nature, and the best check is the self-consistency of the test 
data. 
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Test Results 


Tests were performed with four different ratios of radial clear- 
ance to inner cylinder radius. For each of these configurations 
several series of runs were made with different ratios of outer 
cylinder to inner cylinder angular velocity, and the speed varied 
over the entire possible test range. The temperature difference 
between the cylinders was between 40 and 60 F. 

Fig. 4 and 5 illustrate graphically the flow regimes covered 
for two of the four clearance ratios. On these figures the angular 
velocity of the inner cylinder is plotted against that of the outer 
cylinder. The conditions for each series of runs at a constant 
value of 9 are shown by the dashed straight lines. Superimposed 
on these graphs are the conditions for solid body rotation, a, = 
W:; potential vortex flow, w, = (1 + d/R;)*w.; and the critical 
conditions for the transition from laminar to vortex flow, as given 
by equations (1) and (2). All regions below the latter curve are 
in the laminar flow regime, while above this curve are conditions 
for which Taylor vortexes will exist. Also shown is a curve 
representing a tentatively proposed second transition which the 
present data indicate may exist, and which may separate the 
different vortex regions indicated in Fig. 1. 

The data for the four clearances are shown on Figs. 6, 7, 8, and 
9, where the ratio of Nusselt number to the Nusselt number for 
pure conduction? is plotted versus Taylor number, with cylinder 
velocity ratio as a parameter. For each case, the Nusselt number 
approaches the conduction value, (d/R,)/In(1 + d/R,), for low 
values of Taylor number, and for Taylor number above a criti- 


‘The Nusselt number for pure conduction, (d/Ri)/In(1 + d/R)) 
is a function of clearance ratio only. Its introduction reduces the 
results for different clearances to a common basis for easy comparison. 
The same result could have been obtained by defining the convec- 
tion conductance A in terms of a log-mean heat-transfer area rather 
than the surface area of the heated cylinder. 
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Fig. 10 Correlated test results for inner cylinder only rotating 


cal value the Nusselt number increases. Also shown on these fig- 
ures are the average values of Nusselt number obtained for several 
runs with only the outer cylinder rotating. These data were taken 
principally as a check on the accuracy of the test results, since 
it was not possible, in the speed range attainable, to achieve 
anything beyond laminar flow, or pure conduction behavior, with 
the inner cylinder stationary. 

Fig. 10 shows the data for the inner cylinder only rotating 
(n = 0) for all four clearance ratios tested. Plotted on the basis 
of Taylor number rather than a Reynolds number, all of these 
data fall quite accurately on a single line. Thus it is seen that 
the Taylor number, as suggested by Kaye and Elgar [4], pro- 
vides a good correlating parameter. Above a Taylor number of 
about 90, the data for this case can be well represented by the 
equation 

N Nu 


cond 


Data of Gazley [14] for clearance ratios of 0.0068 and 0.095 are 
also shown on this figure. Probable reasons for the large dis- 
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crepancy between Gazley’s data and that of the present investi- 
gation are presented in a later section of this paper. The results 
of a heat-and-momentum-transfer-analogy solution developed 
later are also included on Fig. 10. 

It would be most desirable to find some combination of varia- 
bles which would bring together the results obtained for the dif- 
ferent cylinder velocity ratios in the same way that Taylor number 
has done this for the different clearance ratios with only the inner 
cylinder rotating. The most obvious requirement of such a 
grouping is that it should bring the transition points, or critical 
Taylor numbers, for all of the cases to one value, or at least to a 
very much reduced range of values. Thus the difference between 
the critical Taylor number at a particular speed ratio, and the 
critical Taylor number for the inner cylinder only rotating might 
be expected to play some part in this grouping. 

Examination of Fig. 6 shows that for this small clearance ratio, 
the horizontal distance between the curves for different speed 
ratios, on the logarithmic scale, decreases greatly with increasing 
Taylor number. The two curves on that figure are, in fact, 
separated by an almost constant difference in Taylor number, in- 
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dicating that for this case it might be possible to reduce the two 
curves to one by using Taylor number minus the difference in 
critical Taylor numbers as the abscissa. 

The data for the larger clearance ratios do not, however, confirm 
this proposal. Examination of Fig. 9, for the largest clearance 
ratio, shows that in this case the horizontal spacing of pairs of 
curves on the logarithmic plot is nearly constant, indicating that 
for this case a promising group to be used as the abscissa might be 
the ratio of Taylor number to the critical Taylor number for each 
of the speed ratios. 

From the foregoing it is apparent that any combination of 
variables which will correlate all of the data for different speed 
ratios as well as different clearance ratios must include the clear- 
ance ratio explicitly, as well as implicitly in the critical Taylor 
number. An empirical combination of variables which reduces 
approximately to the forms proposed for the smallest and largest 
clearance ratios is 


Ny — — 3.5d/R;) 
41.1+(Nn, — R,) 
where N11, is the critical Taylor number at cylinder speed ratio 


n, and Nt, is the critical Taylor number for the inner cylinder 


only rotating, or 7 equal to zero. 41.1 is the value of critical 

Taylor number for 7 equal to zero, and clearance ratio, d/R,, ap- 

proaching zero; this is used rather than Nv... Which is a function 
0 


of d/R,, so that the rather good results obtained by plotting 
against Taylor number alone will not be altered for the case where 
n equals zero. In this grouping, the coefficient of d/R; 
could have any value from about 3 to 4 without making a sig- 
nificant change in the agreement obtained in the results. 

Fig. 11 shows all of the data for concentric cylinders plotted 
on the basis of this group of variables, and falling very close to a 
single curve. In fact, the scatter of the data is little greater than 
that occurring in some of the individual series of runs. For values 
of the abscissa greater than about 2, all of the data may be repre- 
sented to a good approximation by the expression 


— (Naa, — — 3.5d/Ri) 


41.1 + — )(3.5d/R:) 
” erg 


Nx 


“Yeond 


1.1 (5) 
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Final correlation including effects of clearance ratio and rotation of both cylinders in the same or 


The data lie well within a +15 per cent band from this straight 
line. 

The analogy solution is also shown on this graph, since, if the 
proposed combination of variables does actually compensate for 
the rotation of the outer cylinder, the analogy solution should be 
as applicable here as for the case of only the inner cylinder rotat- 
ing. 


Heat-and-Momentum-Transfer-Analogy Solution for 
Concentric Cylinders With the Inner One Rotating 


A heat-and-momentum-transfer-analogy solution for heat 
transfer between concentric cylinders with the inner one rotating 
may be of value both in giving a semitheoretical basis for extend- 
ing the experimental data into other regions, and also in examining 
the flow regimes and transfer mechanisms occurring in this case. 

Gazley [14] has made an analogy solution for very small 
clearance ratios, based only on the assumption of a Karman- 
Nikuradse velocity profile (laminar sublayer and turbulent core 
only) throughout the fluid space, but this assumption actually 
reduces the problem to that of flow between parallel planes in 
motion relative to each other, since it ignores the effects of Taylor 
vortexes in the core, and curvature in general. Gazley’s analysis 
leads to Nusselt number as a function of Reynolds number based 
on inner cylinder velocity and gap width, rather than a function 
of Taylor number, which the experimental results of the present 
work, available friction data, and theories on stability indicate 
to be more pertinent. 

The basic difference between the flow between rotating cylin- 
ders and that between planes in relative motion is the existence, 
even for vanishingly small clearance ratios, of large scale stable 
vortexes in the central part of the annular space between the 
cylinders. In this region the usual concept of small scale turbu- 
lent mixing is not appropriate because of the gross mixing caused 
by the vortexes. Measurements by Pai [6] indicate that the mean 
velocity profile, when averaged with respect to cylinder length, 
has a region of essentially constant velocity in the center of the 
gap, due to the vortex motion. There are also regions near the 
walls in which the velocity distribution is of the logarithmic type 
in agreement with the theory of homologous turbulence. 

In the present analysis, friction data of Wagner [13] and of 
Taylor [12], along with information on transition from laminar 
to vortex flow, were used in conjunction with reasonable assump- 
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tions concerning the transfer mechanisms and velocity profiles 
in the boundary layers in an attempt to predict the heat-transfer 
rate from a heated rotating cylinder to a stationary concentric 
cylindrical shell. 

In essence, the analysis was similar to the usual heat and 
momentum analogies applied to other systems, and the usual di- 
mensionless variables were used, but the form of the velocity pro- 
files dealt with in the boundary layers was somewhat modified to 
fit observed phenomena. The sublayers were considered suf- 
ficiently thin that the shear stress and heat flux through them 
could be treated as constants, and the clearance ratio d/R, was 
considered small enough so that shear stress and heat flux 
through the sublayers on the inner and outer cylindrical surfaces 
could be the same constants. Also, small d/R, allowed the use of 
a simplified shear stress equation, 


du 
= 4 + v) 
Ge dy 


For larger values, du/dy might more properly have been re- 
placed by (du/dy — V/r) for transport of moment of momentum, 
or by (du/dy + V/r) for transport of the forced vortex parameter 
suggested by Kreith [16], the latter form being reasonable since 
this combination appears in the expression for laminar flow viscous 
However, V/r becomes small relative to du/dy with de- 
creasing clearance ratio, and either of the other possible equations 
reduces to the simplified one just given. 


stress. 


The equation for heat flux was 


dt 
(q/A)o = pe,(€y + a) 
dy 
As a basic part of the analogy €, and €y were assumed equal. 
The assumed velocity profile used in the analogy solution 
consisted of a laminar sublayer, a buffer layer, and vortex core. 
The laminar sublayer velocity profile was assumed to be u* = y*, 
and it was taken to extend so far that, at a Taylor number of 41, 
the approximate value of critical Taylor number for small clear- 
ance ratio, the entire fluid space consisted of two laminar sub- 
layers meeting at the center of the space. This, along with the 
laminar flow friction factor, fiam = (2Vd R,) Nv, established 
the limit of the laminar sublayer as y*; = 3.2/(d/R,)'/*. 
A velocity profile of the form used in the Karman-Nikuradse 
profile, but with the constants adjusted to give continuity of slope 
at the junetion with the laminar sublayer was chosen for the 


buffer layer. 
3.2 y* 
‘= 
(a/R) |  3.2/(d/R:) 


Since the measurec velocity profiles of Pai [6], and the whole 
idea of a large seale vortex region with intense mixing, indicate 
that the velocity is constant in the vortex core, the buffer layer 
was assumed to extend far enough so that at its limit, y,*, the 
velocity was equal to half that of the moving inner surface. This 
resulted in the expression 


3.2 (d/R,)'* 
= 7, exp - 


The temperature drops in the sublayers were calculated in the 
usual manner and, since the vortex core has been assumed 
to be a region of uniform velocity and intense mixing, it is thus a 
region of constant temperature. The total temperature drop is 
then only the sum of the drops in the sublayers at the two cylinder 
surfaces. 


Using the friction data on Fig. 2, and setting Prandtl number 
equal to 0.7, leads to the following expressions for Nusselt num- 
ber: 
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For 41 < Ny, < 90 


0.0308 
"1 + 1.43 In [0.7 exp (0.794 — 1) +03] 
For Ny, > 90 
0.0682 
N Nu = (a) 


1 + 1.43 In [0.7 exp (0.358 Ny,°2" — 1) + 0.3) 


These results of the analogy solution have been plotted on Figs. 
10 and 11, along with the experimental results for various clear- 
ance ratios. The curve representing the analogy solution is seen 
to follow the trend of the experimental results quite well, but for 
values of Taylor number above 90 the analogy results are uni- 
formly about 20 per cent above the data points. While closer 
numerical agreement might be desired, the fact that the two 
curves are essentially parallel, and differ as little as they do, lends 
some support to the mechanisms and flow regions postulated in 
the analogy solution. 

A critical review of the assumptions made in the analysis may 
give some clue as to the reason for the discrepancy between 
the predicted and the experimentally observed behavior. 

A fundamental assumption of the analysis was that of the 
equality of the eddy diffusivities for heat and for momentum 
Although this assumption may not be accurate, it cannot be re- 
sponsible for the disagreement, for in the range of Taylor numbers 
up to at least 200, the thermal resistance in the laminar sublayer 
is so predominant that no reasonable adjustment in eddy diffusiv- 
ity could affect the results. 

It is realized that the assumed sublayers imperfectly represent 
the true flow situation especially in the Taylor number range just 
above the transition from laminar Couette flow. In this region 
the laminar sublayer is too thick, since the flow visualization 
photographs of other investigators indicate that only slight), 
above the predicted transition point vortexes develop and occupy 
a large part of the space between the cylinders. At higher Taylor 
numbers the assumed sublayers have a more realistic thickness, 
especially if the vortex pattern proposed by Pai is taken to be the 
one prevailing, and it is realized that the profiles assumed in the 
analysis must attempt to represent an average behavior along 
the length of the cylinders, where vortex motion and small scale 
turbulence seem to alternate in predominance. Thinner sub- 
layers in the low Taylor number range, while perhaps in better 
agreement with observed behavior, would lead only to worse 
disagreement in the predicted heat-transfer performance. 

The final assumption which must be subjected to examination 
is that of a uniform temperature in the vortex region, and it 
appears that in this assumption lies the reason for the discrepancy 
in the results of the analysis. Although the measurements of Pai 
indicate that this is a region of uniform average mean velocity, it 
seems that the vortexes must be an imperfect coupling mechanism 
between the inner and outer sublayers as far as heat transfer i= 
concerned. The fluid in the vortex region might be thought of 
as acting very much like the intermediate fluid in a coupled heat 
exchanger system. If, acting in this capacity, it is not capable of 
transmitting all of the thermal energy which the inner cylinder 
sublayers are capable of delivering to it, then some transfer 
through the vortex region would occur by ordinary conduction 
and eddy diffusion, with an accompanying temperature drop. 
Because of the relative thickness of the vortex region, the tem- 
perature gradient through it would not have to be great to ac- 
count for enough additional temperature drop to bring the results 
of the analysis and the experimental work into correspondence. 
The existence of a small temperature gradient in the vortex 
region is in agreement with a measured temperature profile of 
Taylor [5], which shows a small but continuous temperature drop 
through the central part of the gap between the cylinders. It 
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was thought at first that this observed behavior might be local at 
the particular section of the vortex pattern through which the 
temperature profile was measured, and that it might be compen- 
sated for, on the average, by the temperature distribution at 
another section along the length of the cylinders, as seemed to be 
the case in measured velocity profiles, but perhaps it does repre- 
sent a true average temperature profile. This temperature profile 
was measured at a Taylor number of about 13,000, and with a 
clearance ratio d/R, of about 1.04, both of which are beyond the 
range of the present work, but it does tend to support the reason- 
ing that some temperature gradient may exist in the vortex core. 


Summary and Conclusions 

In this paper test results have been presented for convection 
heat transfer between two concentric rotating cylinders with no 
axial flow of the annulus fluid. The test variables cover clearance 
ratios from 0.054 to 0.246 and six variations of angular velocity 
ratios including counterrotation, corotation, the inner cylinder 
stationary, and the outer cylinder stationary. The flow regimes 
covered are illustrated graphically on Figs. 4 and 5. 

The test results substantiate very well the stability criteria of 
Taylor and Meksyn for transition from laminar to vortex flow. 
It is found that the Taylor number, as suggested by Kaye and 
Elgar, serves very adequately as a correlating parameter for the 
effects of variations in clearance ratio, Fig. 10. It is further 
demonstrated, Fig. 11, that the effects of both clearance ratio and 
angular velocity ratio, including both corotation and counterrota- 
tion, can be correlated by use of a complex parameter based on 
the Taylor number and the critical Taylor number at which 
transition to vortex flow takes place. 

An analogy heat-transfer solution, based on friction measure- 
ments for concentric rotating cylinders, is included to provide 
some insight into the heat-transfer mechanisms involved. It is 
concluded that, at least in the Prandtl number range of air, the 
heat-transfer resistance is predominantly in sublayers close to 
the two cylinder surfaces. However, there is evidence that about 
20 per cent of the temperature change takes place in the vortex 
region, probably due to an imperfect coupling bet ween the Taylor 
vortexes and the sublayers. 
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DISCUSSION 
F. Kreith* 


The careful experimental measurements and the original theo- 
retical treatment of heat transfer through air contained between 
two horizontal concentric rotating cylinders represents a valua- 
ble contribution to the understanding of heat transfer in rotating 
The flow and transfer phenomena in rotating systems 
are exceedingly complex and it appears therefore pertinent to 
point out some possible limitations of the phenomenological 
viewpoint of the authors. Although the following comments bear 
no direct relation to the particular system treated in this paper, 
they may shed some light on transfer phenomena in geometrically 
similar, but fluid-mechanically dissimilar systems. 


systems. 


Since the flow and transfer medium in this study was gaseous, 
it filled the entire gap between the cylinders and was not appre- 
In the analytical treatment of the 
problem the physical properties of the medium were taken into 
account by specifying the appropriate Prandtl number in the 
evaluation of the Nusselt number in Equations (6) and (7) of the 
paper. 


ciably affected by gravity. 


As seen in Fig. 10, the theoretical analogy calculations 
confirm the trend of the experimental data and above the critical 
Taylor number the Nusselt number increases continuously with 
increasing Taylor number. 

In an experimental study carried out at Lehigh University in a 
geometrically similar system, but with water as the transfer 
medium, some additional phenomena were observed. The water 
system differed from the system used by the authors because the 
fluid passed through the annulus at low velocity. The phe- 
nomena observed in the water system were, however, independent 
of the axial velocity and are therefore likely to occur also in an 
enclosed system such as the one used by the authors. 

With both cylinders rotating at the same speed, the Nusselt 
number was found to remain unaffected by rotation until the 
rotational speed reached a critical value. This critical value 
corresponded to a value of the Taylor number given by Equation 
(2) of the paper if the clearance d was taken as one half of the dis- 
tance between the cylinders. Taylor's stability criterion would 
not predict instabilities to occur when 7 is unity, but vortexes 
were visually observed even without any axial flow; i.e., in a 
At Taylor 
numbers above this critical value an increase in rotational speed 
This 


trend continued, however, only until the rotational speed reached 


closed system analogous to that used by the authors. 


produced a corresponding increase in the Nusselt number. 


* Aeronautical Professor of Engineering, University of Colorado, 
Boulder, Colo. Mem. ASME. 
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a value at which the centrifugal force was equal to the gravita- 
tional force. At this point a second break in the curve of Nusselt 
number versus Taylor number occurred. 

Although the reasons for the foregoing phenomena are as yet 
not completely clear, the observations made in the water system 
suggest that in a liquid further complication may be expected 
and the analogy theory may have to be modified. It would 
therefore be of interest to repeat the tests reported in this paper 
using water instead of air as the transfer medium. 

In passing an apparent misprint should be noted. If the forced 
vortex parameter is taken as a transferable property the shearing 
stress will be proportional to (du/dy — V/r), but not to (du/dy + 
V /r) as stated in the paper. 


A.H. Nissan’ and F. C. Haas® 


In Figs. 6to 9 the authors present a number of data which indi- 
It is 
noted that the data represent families of parallel curves of which 
n is the parameter, and the authors proceed to combine the data 
into a single curve as shown in Fig. 11. 

We wonder whether work done by Brewster and Nissan® and 

’ Department of Chemical Engineering, Rensselaer Polytechnic 
Institute, Troy, N. Y. 

¢ D. B. Brewster and A. H. Nissan, ‘The Hydrodynamics of Flow 
Between Horizontal Concentric Cylinders—I, Flow Due to Rotation 
of Cylinder,"’ Chemical Engineering Science, vol. 7, 1958, p. 215. 


cate that the critical Taylor number is dependent upon 7. 


by Brewster, Grosberg, and Nissan’ on instability of flow in 
curved annuli may be of use here. These papers consider the 
onset of instability for a variety of conditions arising from rota- 
tion of cylinders as well as pumping of fluids circumferentially. 
Thus different velocity distributions arise in the gap. In cer- 
tain of these, e.g., when the two cylinders are rotating in opposite 
directions, there will be two distinct regions within the gap sepa- 
rated by a surface of zero velocity. In one of these regions 
circulation will increase outward; in the second, circulation will 
decrease outward. Following Rayleigh's principle, only the 
second region will be potentially unstable; the first is stable. 
Thus instability does not extend over the full gap distance but 
rather occupies only a portion of the gap, designated as 6; 5, the 
width of the region of instability, is easily calculated by an ap- 
proximation for small gaps, 


or more precisely from the Couette flow equation, 


by calculating the radius of zero velocity. Nissan andjhis co- 


7D. B. Brewster, P. Grosberg, and A. H. Nissan, ‘The Stability 
of Viscous Flow Between Horizontal Concentric Cylinders,’’ Pro- 
ceedings, Royal Society (London), in press. 
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Fig. 12 Correlation of Nxu/NNucona with P for all cases given 
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workers define a parameter P which is a modified Taylor num- 
ber based on the instability region, 5, rather than the gap dis- 
tance, and it is defined in a manner which makes it suitable for 
all combinations of rotation and pumping. 


r R, 


P is the average velocity in the gap due to rotation, pumping, or 
both. For the case considered here a simplified expression 


applies. 
(5) Nt 


The factor of one half is introduced to account for the average 
velocity which is required in systems which include pumping. 

On the basis of data taken from the literature for vertical con- 
centric cylinders and their experimental work with horizontal 
concentric cylinders, these workers define the following critical 
values for P: 


0 < b/d 2/; 


/s 1 

d (26/d — 1)? + 0.01 
Applying the parameter, P, to data given in Figs. 6 to 9, it is 
found (Table 1) that P, is indeed fairly constant, having an aver- 
age value of 12.8, for all cases of counter and corotation. The 
difference between the calculated and theoretical values may be 
attributed to experimental error and error in reading the graphs, 
but the constancy of P, is noteworthy. For the two cases of 
corotation in which instability was encountered, the Taylor 
number was adjusted by basing it on the difference in cylinder 
angular velocities, rather than the angular velocity of the inner 
cylinder. The correlation of the given data using P shows itself 
to be quite satisfactory, Fig. 12, herewith. The data lie within 
the +15 per cent band expressed by the authors for the correla- 
tion in their Fig. 11. The data for 7 = 0 similarly fall on one 
line with 


Table 1 Calculated critical values of P 


We 


0.054 
0.084 
0.128 
0. 246 


The one drawback of this approach is the fact of having the 
data fall naturally on two parallel lines. One line correlates the 
data for conditions when both cylinders rotate and the other 
line correlates the data for a stationary outer cylinder; however, 
the correlation is simpler to handle and has theoretical justi- 
fication. 


H. J. Sneck® and G. R. Fox® 


We wish to commend the authors for their excellent piece of 
work and add some additional experimental information which 
may be of interest. 

Fig. 13 shows the results of recent experimental work performed 
at the Bearing and Lubricant Center of the writers’ company. 
These tests were performed on a 4-in-diam bearing having various 
clearance-to-diameter and length-to-radius 
with 70 SUS (at 100 F) oil. The shaft was maintained concen- 
tric with the bearing and only the shaft was allowed to rotate, 
Fig. 14. The oil was fed into the bearing through four feeding 
holes while seals prevented excessive leakage from the ends of the 


ratios, lubricated 


8’ General Engineering Laboratory, General Electric Company, 
Schenectady, N. Y. 

® General Engineering Laboratory, General Electric Company, 
Sehenectady, N.Y. Assoc. Mem. ASME. 
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hearing. Thermocouples were placed in the inlet and exit of 
the bearing and the flow through the bearing was regulated to 
maintain the temperature rise through the bearing less than 80 
The kinematic viscosity 
Total flow 
rates for the oil supplied to the bearing range from approximately 
0.6 gpm for d/h, = 0.004 to about 1.0 gpm for all the other cases. 
The bearing-inlet pressure under these circumstances was in the 


deg F at the highest rotational speeds. 
was evaluated at the arithmetic average temperature 


order of 10 psig. Tn all cases the ratio of the mean axial velocity 
and in most 


eases much less than this value, decreasing with increasing Vr, 


to the maximum tangential velocity was less than 1 
in our experiment. No heat-transfer data were taken, 

The test conditions here differ from those of Wagner and 
Taylor in three important respects 


(a) The medium is a liquid. 


(4) Various degrees of axial-flow velocity are purposely 


introduced 

(c) The length-to-radius ratios are much smaller, correspond- 
ing to the range of values commonly encountered in bearing 
practice 
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For 
a wide range of d/R,, L/R, > 1, and Nr, > 40, agreement is very 
{Nt 


The solid line of Fig. 13 is the same as shown in Fig. 12. 


good. For Nr, < 40, | is independent of Nv, as pre- 


)iam 


dicted by laminar theory. However, when L/R, < 1, the data 


| 


For Nr, > 90 these data exhibit a transition. A closer inspection 
of Wagner's and Taylor’s data as well as the L/R,; > 1.0 data 
of our experiment shows the same tendency. 

Our data suggest that the solid line in Fig. 14 is the asymptotic 
envelope of a family of curves for large L//?; and vanishingly small 
axial to circumferential flow-rate ratio that at higher flow-rate 
ratios these curves extend horizontally to the left (independent 
of Nr,), and that similar families exist for other L/R, ratios with 
the envelope shifting to the left for smaller L//,, Fig. 15. 

We believe that further visual and quantitative investigation 
of these deviations would prove very fruitful. 


show a definite shift upward to higher values of 


Authors’ Closure 


The authors wish to thank the discussers for their interesting 
and valuable comments. 

Although the system described by Mr. Krieth is somewhat 
different from that investigated by the authors, his remarks on 
While the details of his 
tests are not known, it is suggested that because of different ap- 
paratus size, Mr. Krieth’s tests covered an additional flow regime 
which was beyond those which could be obtained with the au- 
thors’ apparatus. 


the phenomena are quite interesting. 


It should be pointed out that the relative im- 
portance of centrifugal and gravitational forces is independent of 
the fluid medium, since both these forees vary in the same man- 
ner with density. 

With regard to the “apparent misprint’? which Mr. Krieth 
points out, it was, unfortunately, not specifically mentioned in 
the paper that the relative velocity u was taken as positive in the 
This defi- 
nition of uv was chosen in order to make the velocity profile equa- 


direction opposite to that of the cylinder rotation. 


tions more similar to the usual pipe flow equations, and with this 
definition the shear stress equations are correct as given. 

The work described by Messrs. Nissan and Haas may indeed 
be of use in providing a correlating parameter for concentric 
cylinder heat-transfer data. It is not obvious to the authors that 
the simplified expression given by Messrs. Nissan and Haas, 


1 /6\"? 
a\a) *™ 


is applicable for any case other than that in which the outer cy!- 
inder is stationary, since this expression implies that 7, the “aver- 
age velocity in the gap due to rotation, pumping, or both,” is 
equal to '/. 


P = 


The authors will be most interested in seeing 
the paper by Brewster, Grosberg, and Nissan when it is available, 
since it will no doubt clarify this matter and also present the theo- 
retical basis for their correlation. 

The additional experimental friction data presented by Messrs. 
Sneck and Fox, and their interpretation of these data, are a valua- 
ble contribution. 
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L. P. HERRINGTON 


The John B. Pierce Foundation, 
New Haven, Conn. 


Full-Scale Human-Body-Model Thermal 
Exchange Compared With Equational 


Condensations of Human 


Calorimetric Data 


INSTRUMENTS Which have advanced knowl- 
edge of biothermal heat exchange are the human calorimeter [4, 
10, 12}! and heated models of the human body [2,5]. Both in- 
struments provide excellent data on a heat exchanger characterized 
by thermal feedback alteration of conductance, evaporative 
process, and heat production. These bio-engineering properties 
combined with complex shape, variable compound insulation, and 
a wide range variation of possible ambient exposures, vastly com- 
plicate the strict application of the rational heat-exchange equa- 
tions established in thermal engineering practice. . 

In an earlier paper [6] prepared for the Committee on Biotech- 
nology of the Heat Transfer Division, the author suggested that 
the empirical human heat-exchange data provided by human 
calorimeters (and models) are most useful to engineers when they 
are condensed into multivariate equations. Work was _ re- 
ported to the committee in which, by algebraic and least-squares 
concepts, a five-element linear differential equation was derived to 
describe 540 hours of calorimetric data on the seated clothed 
human subject in the range of 40 to 80 deg F. The coefficients 
of this equation establish, for this large data collection, the partial 
rates of change of the skin boundary temperature with respect 
to metabolic heat input, evaporative cooling, ambient air, and 
ambient radiation temperature. 

The current report will extend this analysis with respect to four 
principal topies: 


1 The comparison of data obtained with electrically heated 
body models as instruments of measurement with the heat ex- 
change indicated by a linear differential equation for human sub- 
jects in the calorimeter. 

2 The comparison of the original observed calorimeter data on 
human subjects with the heat exchange indicated by the equa- 
tional derivation. 

3. The presentation of four new equations representing the re- 
sponse of different segments of the human body to internal and 
ambient factors affecting heat exchange, and an indication of their 
use in man-machine problems of heat tolerance prediction. 

4 A report of the extension of this study to instrumentally re- 
corded data of the human calorimeter in the region above 80 deg 
F, in which thermally provoked evaporative regulation results in 
an empirically derived linear differential expression of distinctly 
different coefficients. 


' Numbers in brackets designate References at end of paper. 

Contributed by the Heat Transfer Division and presented at the 
Annual Meeting, New York, N. Y., November 30-December 5, 1958, 
of THe American Society or Mecuanicat ENGINEERS. 

Note: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those 
of the Society. Manuscript received at ASME Headquarters, August 
29, 1958. Paper No. 58—A-I181. 
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The Equilibrium Skin-Boundary Temperature of an 
Electrically Heated Model of the Human Body and a 
Vertical Cylinder Compared With Equational Predictions 
Based on Condensed Human Calorimeter Data 

In military supply planning a persistent problem is the rapid 
estimation of the total effective thermal protection of complex 
assemblies of clothing items of changing design and multiple 
combination possibilities. Shape and air-layer factors have made 
the solution of this problem by summation of flat-plate heat- 
transmission resistances unsatisfactory. During World War II 
the use of electrically heated models [2, 12] to determine a re- 
sistance value for each item as worn led to the definition of co! 
units of insulation [3], and the construction of maps showing the 
approximate clo value required at a given season. Thus a map 
location showing a January average minimum of 20 F might be 
indexed as 4 clo. This index value would indicate that a ther- 
mally adequate clothing assembly when placed on an electrically 
heated model would produce a gradient of approximately 64 deg F 
between ambient air and skin boundary as a result of the action 
of a thermal head of about 400 Btu per hr. 

Fig. | shows such a model. Constructed of '/js-in. copper to 
human dimensions with a total surface of 19.40 sq ft (1.8 M?), 
it is provided with segmentally distributed heating units per- 
mitting inputs of | to 4X normal heat production of 400 Btu per 
hr, so distributed as to give segmental skin-boundary tempera- 
tures typical of subjects nude in a comfortable environment at 
85 F, air and wall, relative humidity 50 per cent, and air move- 
ment of 15 ft per min. 

The wide use of such models for tests under a specific condition 
makes it desirable to contrast the values from any generalized 
calculation equation deduced from calorimetric data on human 
subjects with the experimentally observed model values. In the 
earlier paper [6] equation (1) was presented as the best linear dif- 
ferential expression for the relations obtaining between mean skin 
temperature, metabolism, evaporative loss, air temperature, and 
radiant temperature in 180 three-hour calorimeter experiments on 
seated male subjects wearing standard civilian dress (2-piece 
cotton underwear, cotton shirt, socks, low leather shoes, dark 
single-ply suit, */, lined coat, and lined vest). The equation is 
descriptive of the range of ambient temperatures from 40 to 80 F, 
which in physiological parlance is often referred to as the zone of 
body cooling and vasomotor regulation (tissue conductance altera- 
tion by blood flow) [9]. 


NX, = 0.286X. + 0.142NX; + 0.105X, + 0.092X, + 53.39 (1) 
Ts) (T4) (Tw) (M) (E) 


where 7's, T'4, Tw, refer, respectively, to mean skin temperature, 
air and wall temperatures, in deg F; and M and E refer, respec- 
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Fig. 1 Heated model of the human body 


Constructed of '/,,-in. copper and provided with segmentally distributed 
internal electrical heating units which may be adjusted to | to 4 times the 
normal resting metabolism. Copper surface of 19.4 sq ft (1.80 sq m) 
covered with black rubber-base paint and provided with segmental 
thermocouple system registering the mean surface temperature of the 
head, the upper extremities, the trunk, lower extremities, and the total 
skin temperature. 


tively, to metabolism and to evaporation in kg cal per hr per 
man, average height 180 cm (5 ft, 11 in.), average weight 70.8 
kg/(156 Ib), average DuBois surface area 1.76 sq meters (18.9 
sq ft). (In the discussion of the original paper [6], it was shown 
that, in the practical calculation form just given, expression of 
the physiological units in the British system would require 
unnecessary conversion, on each reference to authoritative 
sources from which entry values for M and E must be taken. 
When such tables are common in British units, equation (1) may 
be altered readily to accommodate editorial uniformity. In a 
later section, the coefficients of the equation are given in beta 
form which permits the comparison of the partial rates of 
change in variational form independent of unit systems.) 

We now wish to compare the prediction of the human calorime- 
ter equation (1), with the observations made on the heated model, 
clothed with similar dress, and brought to equilibrium at ambient 
operative temperatures of 60, 65, 70, and 75 F (air and radiant 
temperatures equal at these temperatures and air movement of 15 
ft/min). 
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In Table 1 the observed equilibrium mean skin temperature 
of the clothed copper model has been compared with the mean 
skin temperature of human subjects after 3 hours exposure, as 
computed from equation (1), based on human data from the 
calorimeter. Interpolation in the table shows that at an equal 
air-wall environmental temperature of 68.7 F (15 ft/min) the 
observations on the model and human subjects agree. The cop- 
per model of low thermal] inertia is in complete equilibrium. The 
data on which equation (1) is based involve a mean degree of equi- 
librium for the human subjects after a 3-hr exposure. A reasona- 
ble inference from Table | is that clothed adult males with rest- 
ing heat production about equal to males of average American age 
(cirea 40 years) may be expected to be in approximate thermal 
equilibrium with an operative temperature of 68.7 F after a 3-hr 
exposure. Although heated models are quite convenient for the 
measurement of insulation values, many bio-engineering prob- 
lems of protective clothing require human subjects. The fore- 
going data have value in suggesting an environmental! condition 
and subject status in which human and model evaluation of 
normal clothing may be expected to show agreement in so far as 
the important factor of skin boundary temperature is concerned. 

The thermal theory of static models of such unique shape as the 
human body has not been analyzed thoroughly. However, satis- 
factory ¢xperimental results have been attained in studies of 
radiation area and radiation and air conduction-convection fac- 
tors [7, 8] for the nude human body. 

The experiments reported have afforded an opportunity to com- 
pare in the same laboratory some of the thermal-exchange proper- 
ties of the unclothed model of the human body with the prop- 
erties of a simpler shape. Of the heated shapes that have been 
studied in this laboratory, the one most similar to the approximate 
vertical cylinders which the human extremities and trunk re- 
semble is a blackened vertical cylinder with hemispherical ends, 
24 in. in total height and 8 in. in diameter. The surface emissivity 
of this cylinder is 0.93 per cent. In the calorimeter environment 
in which the model and human data of Table 1 were obtained, the 
heat exchange of this cylinder at surface temperatures 5 to 40 F 
above the environment is accurately given in Btu per sq ft per hr 
by 
H, + H, = 0.51(At)* + 16 X 10- — T,*) (2) 


(H,, H,, convection and radiation exchange in Btu per sq ft 
e r 
per hour; temperature, deg F, absolute scale, second term) 


Table 1 Skin boundary temperatures of blackened copper body model’ 
compared with values from linear differential equation derived from 
calorimeter data on human subjects of comparable size and clothing 
insulation 


(2) 
Observed Predicted 
mean Tain 
equilibrium human (5) 
(1) skin subjects Difference 
Calorimeter temp, (3) 3-hr subject- 
temperature, model, At exposure, model, 
deg F deg F (2 — 1) deg F deg F 
60.1 80.5 +20.4 85.18 +4.68 
65.0 85.3 +20.3 87.28 +1.98 
70.2 90.2 +20.0 89.51 —0. 69 
75.2 95.4 +20 .2 91.65 —3.75 
(68.7)? (88.9) (+20.2) (88.9) (0.0) 


! Surface area of the human body model 1.8 square meters (19.4 sq 
ft), heat input 72.24 kg cal per hour (287 Btu per hour); mean surface 
of the calorimeter subjects 1.76 square meters (18.9 sq ft); predicted 
human subjects skin temperature for the model's heat input computed 
from equation (1), [similar heat input, and insensible evaporative 
loss of 16.46 kg cal per hour per man (65.3 Btu per hour per man) a 
value proportional to the mean relation between metabolism and 
evaporation in the human calorimeter data]. Skin temperature on 
both subjects and model measured at 15 standard points | 12). 

2The environmental temperature at which boundary tempera- 
tures of the model and human subjects are equal, as indicated by in- 
terpolation in Table 1. 
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It is desirable to establish a connection between the heat loss 
of the blackened anatomical model of Table 1 and this more 
common cylinder shape. We report in Table 2 the equilibrium 
surface temperature of the unclothed body model with the same 
heat input (per unit area) as used in the clothed comparison of 
Table 1, and compare the heat input necessary to bring the 
evlinder to the same mean surface temperature. 

The data of Table 2 indicate that, under still air conditions, un- 
complicated by differences in air and radiant temperature, the 
heat loss of a blackened copper cylinder, of dimensions roughly 
comparable to large cylinder-like segments of the human body, 
has combined radiation-convection loss per unit area similar to 
that of our blackened copper model of a human subject. Under 
conditions of forced air movement, different air and radiation 
temperatures, variable clothing insulation, and different postures, 
the comparison is less useful, and many complications have to be 
considered. However, the values of Table 2 do relate the ther- 
mal behavior of the model used in these studies to a heat ex- 
changer of standard shape, under conditions of heat input and 
gradient similar to those which obtain in studies with human 
subjects. 

It is the thesis of this paper and a preceding one [6] that, al- 
though rational calculation methods may be usefully applied to 
analogies existing between cylinder, anatomical model, and the 
nude body in still air, complexities of environment and clothing 
recommend a careful study of equations which condense calori- 
metric data on human exchange with a minimum of assumptions, 


Observational Calorimeter Data on Clothed Human Subjects 
in the Range From 40 to 80 F Compared With a Linear 
Differential Equation Giving the Partial Rates of Change 
of Skin Boundary Temperature as Functions of Metabolism, 
Evaporation, Air, and Radiant Temperature 


Discussion of the data presented by the author at the June, 
1957, meeting of the Society and inquiries relative to the shorter 
form of the paper published in Mechanical Engineering [6] indi- 
cate an interest in the comparison of the predictions of equation 
(1) over the total temperature range of the calorimeter studies. 
This should be presented so that degree of fit of the equation for 
conditions of widely different air and radiant temperature is 
apparent. In the earlier publication [6] the validation of the 
equation as a data fit of 5 variables was indicated by the computa- 
tion of a multiple coefficient of correlation (Ri. = 0.908 + 
0.030), and a comparison of values for the observed and equa- 
tionally predicted means of the five variables studied. 

The intra-range comparisons desired can be obtained from 
Table 3. In this table, 11 sets of calorimeter conditions covering 
180 three-hour exposures on subjects of different build are re- 


ported from the original data. The mean values recorded are the 
clothed averages of subject VII and subject IX for the calorimeter 
conditions indicated, and these data were obtained from Tables 1 
and 2 of the original publication [11]. 

The significant features of Table 3 are column 4, giving the 
mean skin temperature of each group of 3-hr experiments at a 
given calorimeter air and radiant condition; column 5, the pre- 
dicted mean skin temperature from equation (1); and column 6, 
the differences between the observed and predicted skin tem- 
peratures. As may be seen, the validity of the equation as a linear 
differential expression for the partial rates of change of skin- 
boundary temperature with metabolism, evaporation, and ambi- 
ent air and radiant effects, is very satisfactory. For the 11 groups, 
the differences are positive in 6 instances and negative in 5 in- 
stances. The sum of the positive differences for the 11 groups is 
+3.26 F, that of the negative differences —3.23 deg F, vielding 
a mean difference between calorimeter observation and equa- 
tional values of +0.003 deg F. A linear differential expression 
for these data is clearly justified. 

Inspection of the air and radiant temperatures shows that the 
calorimeter data incorporated in the equation cover a wide range 
of air and radiant temperatures, the range of air temperatures 
being from 41 to 72 F, and that of the radiant temperatures from 
41 to 100 F. It is the opinion of the author that all calorimetric 
data on human subjects under various conditions of activity and 
insulation should be condensed in this manner for bio-engineering 
use. It is the intention of the author to accomplish this for all 
data in the Pierce Foundation calorimeter series with which he 
has had an experimental] identification. 

A criticism commonly made of such equations is that, in sub- 
stituting a generalized linear differential expression based on 
physical analysis for specific-case design values and calculation 
from explicit heat-exchange functions, we shorten the process at 
the expense of information on the operation of the separate fac- 
tors. This overlooks the fact that, in animate heat exchangers 
in which properties change on a feedback principle, no efficient 
method of numerical description can easily afford to discard the 
advantages of calculation constants which include the effect of 
the style of response. In short, methods highly reliable for inani- 
mate exchangers of relatively fixed properties are short of efficient 
description if they do not conveniently describe the process of ad- 
justment. It is obvious that equation (1) does include a mathe- 
matical description of the basic physical process and, what is 
equally basic biologically, the process of adjustment as observed. 

In addition to this feature, the process of derivation for such 
equations does permit one to contrast the effective coefficients for 
heat transfer by the various physical processes. This may be 
made clear by the following brief analysis. The coefficients of the 
terms X, — X; in equation (1), its computation form, cannot be 


Table 2 Heat exchange of a blackened copper model of the human body and a vertical hemispherically 
capped cylinder at a common air’ and surface temperature 


Heat 
input, 
Btu 
per hr 
61.1 


Mean 
Area, 
sq ft 
4.19 


deg F 
Vertical cylinder? 89.74 
Human body 

model® 19.4 


287 .0 89.7 


surface, 


Air H,, 
At, temp, Btu/ft? 
Ta — Ts deg F_ perhr 
9.7 80.0 10.06 


H., 
Btu/ft? 
per hr 
4.53 


H, + H., 
Btu/ft? 
per hr 
14.59 

9.7 


80.0 14.7% 


1 Air and wall temperatures equal at 80 F, air velocity 15 ft per min. 


2 Black rubber base paint on ' 


w-in. copper; 24 in. in height, 8 in. in diameter; evacuated and heated 


by an axial Nichrome resistance wire; 8 thermocouples in parallel on outer surface, spaced along the ver- 


tical height. 
5 As described in the text and shown in Fig. 1. 


* Heat input regulated to give surface temperature equivalent to that resulting in the model from 
het input of 287 Btu per hour (72.24 kg cal per hour) per total surface of 1.8 M?. 

5 The partition of the total heat input might be made from radiation area and radiation /convection 
ratios for the human body [7], but direct analysis of these factors for static body models in the standing 


posture have not been reported. 


Journal of Heat Transfer 


aucust 1959 / 189 


gaa 
Prat 
| 
a 
an 
- 
- 
why 
ie 
iy 
sh 
at 
aii 
= 
: 


ter operative temperatures! in the 40 to 80 F range 


Calorime- (4) (5) 
ter (1) (2) (3) Ts T's 
condition To 7 Tw, Exptl, Equat. (1), 
group no deg F deg F deg F deg deg F 
25 $1.72 41.7% 41.7 81.23 80.57 
7 45 5O 45 50.45 83.93 82.93 
9 54.32 32 Si 71.24 81.95 81.50 
Is D4 60 40 19 66 02 81 41 83.69 
5 60 26 59 36 61.07 86.45 8656 
6 61 61 49 46 71.15 84 47 84.78 
7 63 50 41 65 80 60 84 47 S465 
l 7187 71.8% 71.87 91.22 91.16 
5 72 41 49 91 90 14 86 99 86 89 
2 72.41 60.35 81.95 87 98 8833 
4 74.57 41.63 100.49 87.71 867% 


Ta, Te 


, Operative, air, and radiant temperatures, 


directly compared as estimates of the relative weights of the basic 
This is be- 
cause application convenience demands that the variables refer 


processes of convection, radiation, and evaporation. 
both to physical and biological units of measurement. Despite 
the reasonable requirements of editorial order, it is more efficient 
to link the data of fields as separate as engineering and biology, 
retaining the equational entry units which are authoritative in 
the references of the separate fields. Since the engineer applying 
equation (1) will define the environment in deg F, and will use 
tables of heat production and evaporation standardized in kilo- 
gram calories, the equation preserves these units of entry. 

tlowever, we may evaluate the relative importance of the 
variables by reporting here a step in the derivation of the equa- 
tion, 

Equation (1), written in the general form, 


X, = BX, + ByXi+.... (3) 


has coefficients B,, etc., in which both the intrinsic variation of X, 
with the independent variables and the effect of units of measure- 
ment are combined. However, in developing this equation by 
multiple regression methods, deviation units are used, and if we 
set z = (X 
written 


mean)/standard deviation, our equation may be 


= + ... (4) 


By the usual methods these beta coefficients are determined by 
operations such that the average of the squared deviations is a 
minimum. The beta coefficients for equation (1) are reported in 
Table 4 for those who are interested in comparing the relative im- 
portance of the various independent variables in determining the 
skin boundary temperature. 


Table 4 Beta coefficients of equation (1) 


Deviation 


units Relative Process 
82.0 = 0 9395 1 000 Convection 
Buw = 0 7025 0.748 Radiation 
Buws = 0 2492 0 265 Heat production 
= 0.1517 0.172 Evaporation 


The beta coefficients of Table 4 apply to calorimeter heat ex- 
change of the clothed resting subject with mean air temperature 
of 49.03 F + 10.70, mean radiant temperature of 71.52 F + 
16.15, associated with a metabolic heat production of 91.55 kg cal 
per hr + 7.72 (363.3 Btu per hr + 30.6); per man of 1.76 sq m 
(18.9 sq ft); with an evaporative loss of — 20.86 kg cal per hr + 
5.70 (—82.78 Btu per hr + 22.6). 
perature is 61.61 F + 10.17, 


The mean operative tem- 
In general, B. of equation (1) is 
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Table 3 Experimental and equational values for the mean skin boundary temperature of tte clothed human body, male subjects seated at calorime- 


respectively; 


(7) (8) (9) (10) 

(6) Tet, Ts, Metabolism, Evaporation, 

col. 4-col.5 deg F deg F kg cal/hr kg cal/hr 
+0. 66 70 88 98.15 103.52 —16 80 
+100 75.38 98.33 89 60 —15 84 
+0 45 71.33 98.15 100.96 —21 60 
—2 28 71.15 98.15 102.72 —14.7: 
80 78 98.60 88.00 —18.72 
—0 31 77.18 98.51 87.52 —22.24 
-0.18 76.55 98 42 92.48 ~19.52 
+0 06 87 08 98. 64 85.12 —20.96 
4010 98.69 82 40 
—0.35 83.12 98.15 84.80 —31.04 
+0 00 83.03 98 51 90.40 —25 44 

98: 60 —21 


-9 
i“ 


5.70 


Tex, clothing surface temperature; 7's, body temperature (rectal). 


related to 8; (indexed in Table 4 as Bis 45, the partial rate of 
change of X, with NX», 
relation B, = 8,(0,/o2) where sigma is a standard deviation ac- 
companying such mean values as those given previously, the sul:- 
scripts being identified by equation (1). The ratios of the heat 
processes in the experimental data are thus identified. The heat 
input, the temperature gradients, the convection-conduction sur- 
face, and the radiation area for the seated posture are known 
(0.75 X anatomical surface) [7]. These are the elements re- 
quired for a study of the individual processes of heat exchange in 
a physical and analytic sense. It therefore appears to the author 
that this method of condensing calorimetric data for bio-engineer- 
ing use is not limited to the convenience of empirical-data summa- 
tion in a linear differential form but qualifies fully as a technical 
method of advancing the physical description of such complex 
heat-exchange problems by rational methods familiar to the engi- 
neer. 


with variables 3, 4, and 5 constant) by the 


Extension of the Analysis to-the Heat Exchange of the Body 
Segments 

In a previous report [6] to the Committee on Biotechnology of 
the Heat Transfer Division, an equation was given which resem- 
bles equation (1 


except that its purpose is the prediction of the 
local segmental temperature of the legs, as a function of air and 
radiant temperature, metabolism, and evaporation. It was 
thought that such an equation would be useful for those who may 
wish to use an easily measured body-surface temperature as a 
partial index of the state of adjustment of the total surface, or 
in circumstances in which there is a special interest in the effect of 
the environment on a local body area. We now wish to report 
the full set of such equations for the principal body segments 
head, upper extremities, trunk, and lower extremities. 

The various body segments have different properties as heat 
exchangers, and an important element in this is the ratio of surface 
to volume in different body areas. Table 5 gives a good approxi- 
mation of these surface/volume ratios, based upon anatomical 
model measurements. As may be seen, the extremities have a 
large relative surface and are excellent heat exchangers. As such, 
their surface temperature varies much more than the trunk, and 
the environmental condition at which lower extremity surface 
temperatures approach 94 to 95 F is often taken as a topographi- 
cal temperature mark of the inflection of human heat adjustment 
from a phase of adjustment through conductance change to one 
of primary adjustment through sweat secretion [1, 9]. This 
is evident since, as the temperature topography of different 
body segments converges on the 95 F value, additional environ- 
mental heat load is compensated pari passu by evaporation, with- 
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Mean 61 61 05 71.52 85 25 85 25 ais 
std 
dev'n. 10.17 10 70 16 15 3.26 — 
Bes 
4 
- 


Table 5 Surface/volume ratios! of the principal heat-exchanger seg- 
ments of the human body 
Surface, 
per cent 
Body segment of total 
Head 7.0 


Volume, 

per cent Ratio, 

of total surface /volume 
7.0 0 

Upper extremities 18.4 7 72 
Arm 
Forearm 6 10 
Hands 5 3 3.92 


Trunk 31 3 60 


Lower extremities 43 45 
Thigh 18 5 O7 
Leg 16.4 Qf 78 
Feet 7 3.0 2 60 


' The values are considered approximate. Biological variability is 
considerable, and for estimation of surface alone, the DuBois surface 
is generally used, this being based upon a large data collection. In 
the DuBois series the relative surface is head 7 per cent, upper ex- 
tremities 21 per cent, trunk 31 per cent, lower extremities 41 per cent. 
However, associated volume figures were not available. For this 
reason we have used the above comparison based on a body model. 


out marked increases in skin temperature, to the physical and 
physiological limits of the process. 

The surface-temperature response of these various body seg- 
ments for seated, clothed subjects in the range from 40 to 80 F 
operative temperature is known since se smental measurements, as 
well as mean surface temperature, were made in the 180 three- 
hour experiments which were referred to earlier in this paper. In 
Table 6, the observed segmental temperatures for various com- 
hinations of observed air and radiant temperature, metabolism, 
and evaporation are given. In adjoining columns the values ob- 
tained from five element linear differential equations derived to 
condense the data are given. 

A systematic relation exists between the segmental tempera- 
tures and the mean skin temperatures in the calorimeter data re- 
ported in Tables 3 and 6, in that studies of human surface are: 
have indicated the segmental distribution of surface area is as 
follows: Head 7 per cent, upper extremities 21 per cent, trunk 
$1 per cent, and lower extremities 41 per cent [4], and these 
values are established in physiological studies of human heat loss 
as Weighing factors for the derivation of mean skin [7|. Measure- 
ments are generally made at 15 points of definite anatomical loca- 
tion as in the data presentiy being snalyzed 
The following linear differential eq: ations have been deter- 
mined, by methods deseribed earlier in this paper, for each of the 
! primary body segments. In conjunction with equation (1), 
these 5 equations constitute a complete condensation of the 
Pierce Foundation calorimeter studies on seated clothed subjects 
referred to in Table 3. 

The equations comparable to equation (1) for the four principal 
body segments have been analyzed in terms of 7. the operative 
temperature, defined hereafter, and the difference between radi- 
ant temperature and air temperature, in order to study the beta 
co-efficients of air-wall differences, and the relative convenience of 
such expression as compared with separate equational entries for 
air and radiant temperature. The four equations are: 


Ny = O.B15Ns — 0.045X, O.001LX, O.138N 67.00 


Xu 0.550X, — O.144N, + O.085X, OO5LX, 14.96 


Xi. O.138NX, — 0.022N, 0.005X, — 0.006X, 84.55 


= 0.595X, — 0.133X, + 0.206X, + 0.245N 33.61 (8) 


where Xu, Xi, and are, respectively, the temperatures in 
deg F for the head, upper extremities, trunk, and lower extremi 
ties. Ns = 7',, operative temperature [11] technically given by 
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ko 


7. (9) 
where 7'y = radiant wall temperature, 7’, = air temperature, kp 
and ke = subject radiation and convection heat-exchange con- 
stants as determined by calorimeter experiments, and kg = kp + 
ke. Xo = (Tw — and are, respectively, the metabo- 
lism and evaporation in kg cal per hr per man (subjects dimen- 
sionally described in Table 1). The values of variables 8, 9, 4 
are positive, and variable 5 negative, prior to equational entry. 

The coefficients of multiple correlation for the foregoing seg- 
mental skin temperatures are: 7',, head temperature (X), 0.925; 
T,, upper extremities (X),), 0.934; trunk (Xj), 0.880; T,, 
lower extremities (X);), 0.822. These values indicate that such 
linear differential equations fit the original data very well. 
However, a record of the original data and the intra-range fit 
similar to that computed for equation (1) in Table 3 is desirable. 

In Table 6, the mean values of the observed segmental tem- 
peratures for each of the 11 environmental conditions are com- 
pared with the predicted segmental temperatures from the cor- 
responding segmental equations. Inspection of the values of 
Table 6 shows that the intra-range correspondence is good, and 
the mean error small. 


Extension of the Program of Human Calorimeter Data 
Condensation for Bio-Engineering Use to the Environmental 
Region Above 80 F 


In the seated clothed subject at rest or engaged in light ac- 
tivity, the physiological response to increasing environmental 
temperature through the range of comfortable temperatures is 
primarily one of rising skin temperature and increasing peripheral 
tissue conductance, modulated by peripheral blood flow. An 
initial skin-temperature topography under comfortable conditions 
with segmental values ranging from approximately 80 to 93 or 94 
F gradually converges in all segments on a level of 94 to 95 F. 
At this point there is a gradient of 3 to 4 deg F between deep body 
temperature and the skin boundary. Major increase in cardiac 
action would be required to hold deep body temperature constant 
if the process of peripheral evaporation were not stimulated to 
lower radically the rate of rise in skin temperature with increasing 
environmental heat load. 

This temperature control feature of the human body is the fac- 
tor which prevents the extension of equation (1) to air tempera- 
tures of approximately 80 F or higher. 

However, the adjustment in this new range of thermal effect is 
again systematic over a wide range of increasing heat load. In- 
spection of the log of Pierce Foundation human calorimetric 
studies showed that a group of 150 three-hour experiments were 
available ranging in opevative temperature from 79 to 105 F with 
evaporative heat losses ranging from the mean of the value ob- 
served in the sub 80 F series, equation (1), to 125 kg eal per hr per 
man. While 


other data are available, these seemed particularly worth analysis 


This value is some 6 times that of the cooler series. 


since they were based on the same subjects and same clothing 
used in the cool-region analvsis [11 }. 

In the period since the last Biotechnology Committee program 
in San Francisco, these data in the zone of human evaporative 
regulation have been partially processed for the derivation of 5 
variable linear differential equations for mean and segmental 
skin-boundary temperatures. Detailed checking procedures are 
still in progress, and final operations on the correlation matrix for 
However, the means 
check of the observe 1 and equationally pre licted m ean skin tem- 


the segmental data are to be completed. 


Since it is a reference convenience 
to have comparable equations for the 40 to 80 F and the 80 to 105 
F regions in the same source, the new equation is reported here 


perature is quite satisfactory. 
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Table 6 Experimental and equetional values for skin —— temperatures of the 4 principal body seg- 


ments. Group for |, male 


40 to 80 F range. 


Calorime- 


bjects, at calorimeter operative temperatures! in the 


(5) (6) (8) 
ter (1) (2) ( 3) (4) Exptl. Equat. (5) & Exptl. 
condition To, Ta, Tw, od Ta, Tu, TH. Jol. Tu, 
groupno. deg F deg F deg F deg F deg F deg F 5-6 deg F 
25 41.72 41.72 41.72 0.00 81.59 84.16 —2.57 77.27 
17 50.45 50.45 50.45 0.00 86.09 85.08 +1.01 80.87 
19 54.32 32.81 71.24 38.43 85.01 85.36 —0.35 77.09 
18 54.60 40.19 66.02 25.83 84.29 85.07 —0.78 76.64 
5 60.26 59.36 61.07 1.71 89.24 88.48 +0.76 83.57 
6 61.61 49.46 71.15 21.69 88.79 88 50 +0.29 80.51 
7 63.50 41.63 80.60 38.97 88.70 87.94 +0.76 80.60 
1 71.87 71.87 71.87 0.00 92.05 92.52 —0.47 91.05 
3 72.41 49.91 90.14 40.23 92.03 91.34 +0.69 85.10 
2 72.41 60.35 81.95 21.60 91.67 93.12 —1.45 87.17 
4 74.57 41.63 100.49 58 .86 91.22 91.35 —0.13 85.55 
Mean 61.61 49.03 71.52 22.48 88.24 88.45 82.31 
Std. de- 
viation 10.17 10.70 16.15 19.39 >. 67 4.54 
Difference, 
Eexperimental-equational value —0.21 
Calorime- (9) (11) (12 14 (15) 
ter Equat. (6) Exptl (7) (13) Exptl. Equat.(8) (16) 
condition Tv, ol. Tr, Col. Col. 
group no. deg F 8-9 deg F deg F 11-12 deg F deg F 14-15 
25 75.84 +1.43 90.32 89.89 +0 43 76.28 75.64 +0. 64 
17 79.52 +1.35 91 31 91.16 +0.15 79 61 78.21 +1.40 
1” 76.78 +0.31 80 87 90.83 —0.96 77.90 76.63 +1.27 
IS 79.25 —2.61 90.50 91.10 —0.60 76.55 80.21 —3.66 
5 84.38 —0.81 92.84 92.50 +0. 34 82.85 82.77 +0.08 
6 82.04 —1.53 92.30 92.26 +0.04 79.88 79.97 —0.09 
7 81.04 —0.44 92 66 92.11 +055 79.61 80.48 —0.87 
l 00. 56 +0.49 O01 94.17 —0 16 S807 88.83 +0.14 
84.77 +0 33 93.38 93.39 —0.01 82.31 82.39 —0 O8 
2 87.31 14 93 11 93 83 —0.72 83.68 +0.43 
83.87 +1.68 93.83 93.25 +0. 5S 83.66 82.53 +1.13 
Mean §2.31 92.19 92.23 81.07 81.03 
Std. de- 
viation 1.49 4.35 
Difference, 
Eexperimental- —0.00 —0 04 +0.04 


equational value 


1 To, Ta, Tw, operative, air, and radiant temperatures, respectively; 
upper extremities, trunk, and lower extremities, respectively 
ration values needed in conjunction with columns 2 and 3, for calculation of columns 6, 9, and 


12, are given in Table 3 


Table 7 Mean values of hot-region calorimetric data, mature male subjects, 


Tu, Tv, Tr, Tx, head, 


Additional metabolism and evapo- 


normal clothing, seated 


posture 
Skin Air Radiant Operative 
temp, temp, temp, temp Metabolism, Evaporation, 
deg F deg F deg F deg F kg cal/hr kg cal/hr 
Mean 94.23 75.60 105.73 92 43 83.88 —67 .02 
Std. deviation 2.00 12 67 12.75 7 7.10 34.63 
No. of experi- 
ments 150 150 150 150 150 150 


' Subject characteristic 
this report. 


The associated regional analysis and detailed data comparison, 
comparable to Tables 3 and 6, are deferred until this study is ex- 
tended to the analysis of the humidity and vapor-pressure fea- 
tures of similar data. 

In Table 7 are the general descriptive features of equation (10) 
now available as a condensation of human calorimetric data on 
the seated clothed subject for the range from 80 to 105 F. 

The linear differential description of the data of Table 7 is 
given in equation (10). 


= 0.171. X + 0.137X + 0.082Xi7 + 0.02X + 61.34 (10) 
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8, units, and variables identical with those reported in [6] and other tables of 


Equation (10) gives a value for mean skin temperature of 94.29 
F when the means of the independent variables are used to enter 
the equation. The experimental mean of the series is 94.23 F. 
This is satisfactory as a means check, and will be supplemented 
in later work with an intra-range comparison. Since the multiple 
correlation coefficient of equation (10) is +0.93 it is probable that 
the intra-range fit will be of the same order as that reported in 
Tables 3 and 6 for the region below 80 F. 


Summary 


It has been the purpose of this paper, considered jointly with 
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a preceding one [6], to serve several reasonable interests of the 
Committee on Biotechnology of the Heat Transfer Division. 


(a) One of these aims has been achieved by condensing calori- 
metric data on seated, clothed human subjects into linear differ- 
ential equations which permit the bio-engineer to quickly ap- 
proximate the effect of complex air and radiant temperatures on 
the human element of a man-machine design problem. 

(b) These aims have been accomplished for a range of thermal 
loads typical of the region from 40 to 80 F, and separately by 
reason of different human temperature regulation characteristics, 
for the range from 80 to 105 F. 

(c) Static heat exchangers—an electrically heated model of 
the human body and a hemispherically capped eylinder—have 
been used to establish an experimental connection between the 
heat-exchange properties of these inanimate bodies and linear 
differential equations which generalize the thermal interchanges 
of the human body as observed in a calorimeter environment. 
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DISCUSSION 
E. F. Adolph? 


The author reduces the results of certain measurements of heat 
exchanges by clothed men to linear equations. These equations 
are based on the concept that the physical components of ex- 
change—radiation, conduction, convection, and evaporation—as 
well as metabolic heat production, combine to influence the mean 
surface temperature of the human body. Radiative or wall tem- 
peratures were varied independently of the convective or air 
temperatures, to provide a wide range of atmospheric conditions 
for men from whom the necessary data were obtained. A previous 
paper of the author covered the range of ‘operative’ environ- 
mental temperatures from 40 to 80 F; the range 80 to 105 F is 
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now included by a similar linear equation which, however, has 
different coefficients for each of the five terms. In this range the 
physiological regulation of body surface temperature differs 
greatly, since evaporative cooling now predominates. 

It is very satisfying to have an empirical equation that con- 
denses into it the varieties of information needed for assessing the 
steady state of heat regulation achieved by resting clothed in- 
dividuals. The author’s contribution consists chiefly in the use of 
surface temperature as the parameter most readily related to the 
physical conditions of environment. By this means, his two equa- 
tions, of similar dimensions, become available to describe the net 
physiological! state in any environment within the wide range in- 
dicated. 

I agree with the author that this method of presenting calori- 
metric data advances the physical description of heat exchanges 
of the human body. It would be interesting to compare the ac- 
curacies of prediction of physiological effects, first from these two 
equations, and then from charts for effective temperatures and 
for other empirical relations that have been in use. 

I also like the author’s partial evaluation of the differences in 
heat exchange between the “copper man’’ as a model system and 
the average human subject. Differences between model and 
subject will allow us to identify special characteristics of heat 
regulation in man. 


Alan H. Woodcock® 


This paper, extending an earlier one by the same author, will be 
welcomed both by engineers and others who must consider the 
thermal relationships between man and his environment. The 
author and his co-workers have accumulated great quantities of 
pertinent information. However, to those not thoroughly familiar 
with this field, such masses of data are often baffling. His re- 
duction or condensation of it into the form of simple linear equa- 
tions, which can be easily used, is indeed a real service and, need- 
less to say, must be performed by one thoroughly familiar with it. 

The results with the life-sized model of man in this paper are 
reassuring as they indicate that man obeys the physical laws of 
heat transfer in a manner similar to an inanimate object. Table 1 
indicates that the equation derived from human data gives good 
agreement with the model at an air temperature of 68.7 F. How- 
ever, over the relatively small range (15 deg) of calorimeter tem- 
perature, the difference in skin temperature between equation 
and model shows a rather poor correlation. For the model with 
fixed heat input, the term in equation (1) representing metabolic 
heat input is constant and that for evaporation is also considered 
asaconstant. Reference to Table 3 shows that neither metabo- 
lism nor evaporation remain constant for man as ambient tem- 
perature varies. If it did, these terms would have been omitted. 
Had the author used a variable heat to his model to correspond 
with man’s variable metabolic heat production at varying tem- 
peratures and also corrected his evaporative term to correspond 
with temperature, I am sure he would have obtained a closer agree- 
ment between model and equation at al] temperatures. Such 
corrections are, however, complications which are not necessary if, 
as the author recommends, the correct ambient temperature is 
used. 

His introduction of a much simplified model in the form of a 
cylinder with hemispherical ends has both advantages and dis- 
advantages. Its advantages are the reduced size and simplifica- 
tion of shape which, of course, reduces cost and at the same time 
makes an instrument much easier to store and use. Its disad- 
vantage is that the effects of different types of clothing cannot be 
compared on the cylinder since it is impossible to duplicate the 
fit and drape of clothing. However, such complications are 
really beyond the scope of the author’s objective which is to rid 
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the problem of its complications and present the data in its original 
form. An instrument somewhat similar to this cylinder, which 
is used in a slightly different manner, is the globe thermometer of 
Vernon. This globe thermometer, which is a six-inch diameter 
hollow-blackened copper sphere with a thermometer inserted in 
it, has been used with considerable success to integrate the effects 
of air and radiant temperature and is rapidly gaining popularity 
in this country. 

In a later section of his paper, the author proposes a new equa- 
tion [his equation (10)] to extend the range up to 105 F. A com- 
parison of equations (1) and (10) shows that all the numerical 
coefficients in equation (10) except the constant term have bee. 
reduced. This indicates that at higher temperatures variation 
Indeed, skin temperatures be- 
come a relatively poor index of heat stress on man since they are 


of skin temperature decreases. 


so Insensitive 

The effect of humidity is also not included in equation (10) al- 
though it is well known that it is an important factor. However, 
he indicates that his work will be extended to include these fae- 
tors. Such an extension is definitely required and it is hoped that 
he will complete it in the near future. 


Author's Closure 


Dr. Woodcock comments on an interesting feature of equation 
(1) for the 40 to 80-deg F range as compared with equation (10) 
for the 80 to 105-deg F range. He notes that the coefficients of 
equation (10) for the hot region are reduced with respect to equa- 
tion (1) for the cool region. This fact is interpreted as an indi- 
cation of a smaller skin temperature variability at higher ambient 
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temperatures, and possibly as a restriction on the use of skin tem- 
perature as an index of heat stress in the hot region. It is true 
that the rate of change of skin temperature is slower as ambient 
stress increases. However, where the measurement panel con- 
tains the evaporative value whose systematic increment (of op- 
posed sign) maintains this slower skin rise, there is no change in 
the total amount of information incorporated by empirical mul- 
tiple variable equations of equally good fit in the two regions. 
In general, one cannot compare the absolute coefficient values of 
such equations as (1) and (10) without expressing them in beta 
form, see equation (3). 

Likewise, the equations do not contain a thesis as to which 
variable is the best index of heat stress. From the value of the 
multiple correlation coefficient for equation (10) the total infor- 
mation extracted from the hot-series data is of the same order as 
that incorporated in the cool-region analysis equation. 

These equations are numerical devices for expressing condensed 
interrelations as measured. Since they have simple algebraic 
properties, the total evaporation (£) could be made the depend- 
ent variable, if we regard this as the best thermal stress measure. 
In bio-engineering application within the range of the equation 
and its individual variables, the data interrelations will remain 
constant for this body of calorimetric data, analyzed hy the meth- 
ods described, without regard to which variable is regarded as the 
predicted variable. 

An analysis in process will include a vapor pressure term. For 
documentation of equation (10) which is briefly mentioned in the 
text, and equation (1) it may be added that relative humidities 
were held between 40 and 50 per cent for every experimental air 
temperature. 
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The Emissivity and Absorptivity 
of Parachute Fabrics 


The use of parachutes for recovery of information and equipment from high-speed 
vehicles has directed attention to the heating problem which in some instances may be 
so severe as to cause destruction of the parachute. Consequently, the parachute engineer 
requires sufficient information on the heat-transfer characteristics of geometries re- 
sembling those of parachutes and on the heat-transfer properties of the materials used to 


allow a calculation of the temperature history of a descending parachute. In particular, 
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the radiation properties of the parachute materials must be known since the parachute is 
receiving radiant energy from the sun and loses energy by radiation to the surroundings. 
The measurement of radiation properties for such parachute materials is more compli- 
cated than for solid surfaces since an appreciable amount of energy is transmitted by 


such fabrics. It is the purpose of this paper to describe the equipment which was used 
to measure the absorptivity for sclar radiation and the long-wave-length emissivity for 
such parachute materials and to report these data for a number of important parachute 


materials. 


Riaceue a material in the form of a thin sheet as 
shown in Fig. 1 and radiant energy consisting of parallel or nearly 
parallel rays impinging on surface 1 of this material under an 
ungle @. Part of this radiation is reflected from the surface. 
Two limiting cases for the space distribution of the reflected radia- 
tion are specular and diffuse reflectivity. Usually, the actual re- 
flection of surfaces used in engineering is in its space distribution 
between these extremes with larger amounts of energy reflected 
in a range near the specularly reflected rays, but some reflection 
occurring diffusely also [1].? In its details, the space distribution 
depends on the surface condition and on the material of the sheet. 
In this paper, the space distribution of reflected radiation will not 
be investigated. The ratio of the total amount of reflected 
energy—regardless of direction—to the amount of impinging 
energy will be called reflectivity and will be denoted by p. 

The radiant energy which penetrates through the surface 1 is 
partially absorbed in the sheet. The rest leaves it through its 


surface 2. The directional distribution of this radiation depends 


' Publication from the Heat Transfer Laboratory, University of 
Minnesota, Minneapolis, Minn. 

2 Numbers in brackets designate References at end of paper. 

Contributed by the Heat Transfer Division and presented at the 
Annual Meeting, New York, N. Y., November 30-December 5, 1958, 
of THe AMERICAN Society OF MECHANICAL ENGINEERS 

Nore: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those of 
the Society. Manuscript received at ASME Headquarters, July 31, 
1958. Paper No. 58-—A-125. 


Nomenclature 


again on the surface condition and on the material of the sheet. 
The ratio of the total energy integrated over all directions in 
space—which leaves surface 2 to the impinging energy will be 
called transmissivity and will be denoted by r. 

The rest of the radiant energy is absorbed in the material so that 
with an absorptivity @ defined in an analogous manner, the follow- 
ing equation holds 


pta+re=l (1) 


The radiation properties a, p, and 7 depend generally on the 
temperature of the sheet and on the wave length of the impinging 
radiation. In engineering applications, the impinging radiation 
comprises usually a large range of wave lengths and the radiation 
properties consider the respective energies integrated over all 
wave lengths. It has to be kept in mind that these total radiation 
properties depend then on the wave-length distribution of the in- 
cident radiation. In particular, energy considerations on para- 
chites have to deal with solar radiation and with radiation from 
sources of moderate temperature. The radiation properties in 
this paper consider, therefore, these radiation sources. 

According to Kirchhoff’s law, the monochromatic emissivity 
is equal to the monochromatic absorptivity. For the values which 
are integrated over all wave lengths, however, the relation 


awe 


holds only when the source for the incident radiation is a black 
body with a temperature equal to the temperature of the absorb- 


amount of energy reaching = 
thermopile from. test sur- T 
face divided by actual 
amount leaving test surface @ 
amount of energy emitted by 
surface, Btu/hr sq ft 
radiation impinging on ther- 
mopile receiving 
Btu/hr sq ft 
incident energy, solar energy 
or black-body radiation, 
Btu/hr sq ft 


total 


surface, 
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proportionality constants 

= absolute temperature, deg R 
total absorptivity 
total emissivity 

hemispherical 
tivity 

reflectivity for radiation emit- 
ted by thermopile surface 


total hemispherical transmis- 
sivity 


= Stefan-Boltzmann constant, 
0.1714 Btu/hr ft? 
A = galvanometer reading 


Subscripts 
reflec- = hemispherical mirrored surface 
= refers to readings when reference 
surface is in position 
refers to readings when test surface 
is in position 
thermopile 
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Fig. 1 Radiation reflected and transmitted on a thin sheet 


ing sheet, and as good approximation, where the two tempera- 
tures do not differ much.* In this paper, the absorptivities pre- 
sented in Table | can be interpreted as emissivities. 

It can be expected that the radiation properties as just defined 
change with the angle of the incident radiation. This has been 
observed for solid materials of a thickness sufficient to make the 
transmitted radiation practically zero [2]. For thin sheets like 
the materials investigated in this paper, an even stronger de- 
pendence on angle may be expected. The measurements re- 
ported here are all for normal incidence of the radiation. 

The radiation properties depend also on the medium which sur- 
rounds the sheet. In the present experiments this medium was 
always air. 


Description of Apparatus 


The hemispherical integrating radiometer, shown in Fig. 2 and 
previously described in [3], the reflectivity 
and transmissivity of the various parachute fabrics for normal in- 
cident solar radiation. 


was used to measure 


These radiation properties were determined 
also for infrared radiation using the same basic instrument, but 
with the source of radiant energy being a low-temperature black 
body as shown in Fig. 3. For convenience in the program, a 
separate radiometer with the minor modifications shown in this 
figure was constructed for the infrared measurements. For both 
instruments, the hemispherical shells are made out of glass with 
9.625 in. ID and with apertures 9/16 in. diam. The inner surface 
of the hemisphere is coated with aluminum selected because of 
its resistance to tarnishing and also for its relatively constant 
and high spectral reflectivity. The hemisphere is surrounded by 
a cooling jacket to maintain defined conditions around the ther- 
mopile and test surfaces at all times. 

A Kipp and Zonen Solarimeter, type G-19, thermopile without 
base and screen is used in the instrument with the receiving sur- 
of the thermopile The 
and the mounted on a stand which 
can be rotated so that radiant energy is directed either on the 
thermopile or on the test surfaces. 


face being approximately 1/4 sq in. 


thermopile test surfaces are 


The test-surface holder for the 
solar measurements is provided with positions for six test sur- 
faces. 


When the measurements on one surface are completed, 


another surface can be rotated into the measuring position. 


* An estimate of the deviation of both values can be made from 


Fig. 108 of reference [2}. 


196 


‘ AUGUST 1959 


Table 1 Radiation properties of parachute cloth for normally incident 
black-body radiation at 350 F 

Reflec- Trans- Absorp- 
Material tivity missivity  tivity* 
Dacron, 100 lb... 0.08 0.17 0.75 
Dacron, 300 lb.... 0.07 0.14 0.79 
Dacron, 600 Ib... ... 0.07 0.07 0.86 
Dacron, 800 lb.. 0.07 0.03 0.90 
Nylon rip-stop 
(orange) 1.1 oz per sq yd, t- 0.05 0.19 0.76 
MIL—C—7020 } 
Nylon rip-stop ) 
1.1 oz per sq yd (white) } 0 06 0.21 0.73 
MIL—C—7020 
Nylon rip-stop ) 
1.6 oz per sq yd (white) >. 0.07 0.21 0.72 
MIL—C -7020B Type Ill} 
Nylon cloth 2.25 oz per sq yd, 
MIL—C—7350B Type I.. 0.07 0.05 0.88 
Nylon cloth 4.30 oz per sq yd, 
MIL—-C—-8021 Type I 0.07 0.06 0.87 
Nylon cloth 7.0 oz per sq yd, 
MIL—C—8021 Type I. 0.10 0.05 0.85 
Nylon cloth 14.0 oz per sq yd, 
MIL—C—8021 Type 1. . 0.10 0.04 0 86 


® This value is equal to the emissivity of the material at 350 F 
when the monochromatic absorptivity is considered independent of 
temperature. 


Since the fabric test surfaces transmit energy, the background 
behind the test surfaces for the reflectivity measurements was 
designed to absorb any transmitted energy. For the solar-energy 
measurements the background surface was the blackened top of 
the test-surface holder, while for the infrared measurements, the 
background surface was a water-cooled black body as shown in 
Fig. 3. The energy reaching the thermopile when a test surface 
was in test position over that coming from the background was , 
then due to the reflected energy. 

When making the solar measurements, a fused-quartz window 
2 mm thick was installed over the aperture to eliminate effects of 
convection currents on the thermopile readings. Throughout the 
important part of the solar spectrum, 0.3 to 3.0 yw, the fused 
quartz filters evenly and has a transmissivity of 92 per cent [4]. 


Experimental Procedure 


A comparison technique is used for the measurement of the re- 
flectivity. The thermopile and the test surface are located on 
conjugate foci of the hemisphere. If radiant energy is reflected 
then the necessary integration of energy re- 
flected at all angles is accomplished by the hemispherical surface 
since, 


from the test surface, 


for the small distance (0.75 in.) between the foci, the hemi- 
sphere is practically identical in shape with an ellipsoid. The read- 
ing of the galvanometer connected to the thermopile is then an in- 
dication of the amount of energy reflected at all angles since it 
is integrated by the hemispherical surface and focused onto the 
thermopile. By taking the ratio of the energy reflected from the 
test surface onto the thermopile to that reflected by the reference 
surface on the thermopile and using the value of the reflectivity 
of the reference surface, the reflectivity of the test surface is cal- 
culated. The details of the data analysis may be found in the 
Appendix. 

To obtain a high accuracy with the described technique, the 
reference surface should have reflection characteristics similar to 
the test surfaces; i.e., similar directional distribution of the re- 
flectivity. In the visible region of the spectrum one can deter- 
mine readily whether the reference surface and test surface reflect 
similarly, while in the infrared one should use a reference surface 
made from the same type of material. A reference surface which 
satisfies these requirements in the solar region of the spectrum is 
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Fig. 2 Hemispherical integrati 
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THERMOPILE LEADS 


magnesium oxide (MgO). The MgO surface was prepared by 
methods for which the reflectivity was well established [5-10]. 
The reference surface consisted of an MgO layer, 1 to 2 mm thick, 
smoked on a block of compressed MgO. This surface has a re- 
flectivity in the solar spectrum of 0.94 + 0.02. The reference 
surface used in the infrared region of the spectrum was a cloth 
material designated QM7 (Quarter Master Corps standard cloth) 
for which monochromatic reflectivity information was available 
in the literature [11]. A check on the reliability of the infrared 
reflectivity measurements was made by measuring the reflectivity 
for normally impinging solar radiation and by comparing it to 
the value calculated from the information in this reference. The 
value determined in this investigation was 0.29, while the value 
obtained from [11] was 0.25. The reflectivity of the reference 
surface for black-body radiation from a source at 350 F was then 
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for solar-energy measurements 


calculated from the information in [11] and found to be 0.12. 

For the transmissivity measurements, the test surface holders 
were modified so that the test surface, when in the measuring 
position, completely covered the receiving surface of the thermo- 
pile. With no test surface in position, the radiant energy passing 
through the aperture was made incident on the thermopile. The 
thermopile indication is then a measure of the incident radiant 
energy. 
Since the test surface is directly over the thermopile, all the 
transmitted energy was reaching the thermopile. By taking 
the ratio of the transmitted energy to that directly incident 
on the thermopile, the transmissivity is determined. 

From these measured values of the reflectivity and transmis- 
sivity, the absorptivity of the test surface is calculated with 
Equation (1). 


Next, the test surface was rotated into the test position. 
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Fig. 3. Hemispherical integrating radiometer for infrared radiati t 


Test Results and Discussion 


Radiation Properties for Incident Solar Radiati To establish the 
reliability of the present measuring technique, data were first 
obtained using cloth samples supplied by Profs. R. V. Dunkle 
and J. T For these ma- 
terials, the radiation properties are reported in [11]. 


Gier of the University of California. 
Table 2 com- 
pares the results of our measurements with those in this reference. 
It may be seen that they are in substantial agreement. 

Having established confidence in the measuring technique, the 
measurements were then extended to the parachute fabrics 
The resulting values of the reflectivity and the transmissivity for 
solar radiation are given in Table 3. 
lated from Equation (1). 


The absorptivity is ealeu- 


tepresentative photomicrographs of the parachute fabrics are 


shown in Fig. 4. To allow for possible variations in the materials, 


Table 2 Comparison of present reflectivity values with previous in- 
vestigations 

Measured by 
reference {11 | 


Present 
Material measurements 
QM1, cotton sheeting 
bleached, 4 oz per yd 
QM2, cotton sateen 
prepared for dyeing, 9 0z pervd 0.68-0.72 0.70 
QM4, cotton sateen, 


0 62). 66 0.63 


undyed, 9 oz per vd 0. 69-4). 72 0.69 
QM6, cotton sateen, 

medium gray, 9 oz per vd 0.53 0.54 
QM7, cotton sateen, 

dark gray, 9 oz per yd 0.24 0.29 
50°, wool, 50°, cotton, 

knit, undyed, 10.5 oz per yd 0.62 0.60 
Cotton knit, undyed, 

3 oz/vd 0.60 0.58 
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Table 3. Radiation properties of parachute cloth for normally incident 
solar radiation 


teflec- Trans- Absorp- 
Material tivity missivity tivity 
Dacron, 100 Ib. 0.35 0.60 0.05 
Dacron, 300 Ib 0 54 0.35 0.11 
Dacron, 600 Ib 0 61 0.27 0.12 
Daeron, 800 Ib 0 62 0.19 0.19 
Nylon rip stop 
(orange) 1.1 oz ner sq vd. > 0. 23 0 64 0.13 
MIL—-C—7020B Type I 
Nylon rip-stop 
1.1 oz per sq vd (white) > 0.27 0.65 0.08 
MI'.-C—7020 } 
Nylon rio-stop 
1.6 92 per sq vd (white ie 0.22 0.72 0.06 
Type III 
Nylon cloth 2.25 oz per sq vd, 
MIL--C——-7350B Type I 0.36 0.59 0.05 
Nylon cloth 4.30 0z per sq yd, 
MIL—C-—-8021 Type lL. . 0.44 0 48 0.08 
Nvlon cloth 7.0 0z per sq vd, 
MIt--C—-8021 Type I 0.46 0.41 0.13 
Nylon cloth 14.0 0z per sq yd, 
MIL—C— 8021 Type III... 0 #2 0.27 0.11 


several different samples were used for each material and the re- 
sults were averaged. It is apparent from Table 3 that the major 
portion of the impinging solar radiation is either reflected or trans- 
mitted, and that the resulting absorptivity values vary from 0.05 
for the lightweight weaves to 0.19 for the heavier fabrics. 
Radiation Properties in Infrared Region of Spectrum. The measuring 
technique used in this study was analogous to the one for the 
solar-radiation measurements. The hemispherical radiometer 
with the black-body source is shown in Fig. 3. Using the ref- 
erence surface QM7 and the test surface, the reflectivity and 


Transactions of the ASME 


N 
THERMOPILE 
(Aad 
| AX 
N 
ALUMINUM 


DACRON 


Fig. 4 Parachute fabrics micrographs. Dacron—(c) CO Ib tensile strength, (b) 600 Ib tensile strength, (c) 
300 Ib tensile strength, (d) 100 Ib tensile strength. Nylon cloth—(e) 14 oz per sq yd, MIL-C-8021 type Ill; 
(f) 7 oz per sq yd, MiL-C-8021 type Il; (g) 4.30 oz per sq yd, MIL-C-8021 type |; (h) 2.25 oz per sq yd, MiL- 


C-7350B type I. 


Nylon rip stop—(i) (orange) 1.1 oz per sq yd, MIL-C-7020B type |; (j) (white) 1.1 oz per 


sq yd, MIL-C-7020; (k) (white) 1.6 oz per sq yd, MIL-C-7020B type Ill. 


transmissivity measurements were made. Again the absorptivity 
was determined from these two values and the results are col- 
lected in Table 1. 

A black body at a temperature of 350 F was used in this series 
of tests, and consequently, the impinging radiation was concen- 
trated in the far infrared portion of the spectrum. As contrasted 
with the solar portion of the spectrum where the absorptivity 
values are low, the infrared absorptivities are appreciable ranging 
from 0.72 to 0.90. 

If the monochromatic absorptivity is assumed to be invariant 
with temperature then this measured absorptivity is equal to the 
emissivity of the parachute material at 350 F. 
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APPENDIX 
Analysis of Data 


Reflectivity. The radiant energy G arriving at the thermopile 
! when the test surface s is in position with incident radiation 
J, impinging on the test surface may be written (Fig. 5) 


The first term on the right represents the incident solar or infrared 
energy reflected off the test surface corrected for the reflectivity 
of the mirrored surface and the loss through the aperture; the 
second term represents the energy emitted by the test surface 
and its background with the same corrections, while the third 
199 
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term is the energy originating at the thermopile itself which ar- 
rives back on the thermopile after reflection from the test surface. 
If the aperture is now closed to external radiation so that J, is 
zero, the energy arriving at the thermopile becomes 
G(2s) = e,p,,C, + eC ,7p,,7p,' (3) 
Solving Equations (2) and (3) for the reflectivity of the test sur- 
face p,, 
— G(2s) 


J Pf 


(4) 


The same type of relation obviously holds for the reference surface 
(Tr) 
Gir) — G(2r) 
= 


(5) 


r 


By a combination of equations (4) and (5), the reflectivity p, is 


obtained 
C, \ GAs) — 
( G(r) — G(2r) 

A heat balance may now be written on the thermopile itself 
when it is receiving energy in the amount G(1s), since at equi- 
librium this energy must be balanced by the energy lost by the 
thermopile surface by radiation, conduction, and convection: 


(6) 


= + ki(Ty — T.) (7) 
where the first term on the right represents the radiation loss and 
the second term the convection and conduction losses. These 
losses can be assumed proportional to the difference between the 
temperature of the receiving surface and the temperature of the 
cold junction of the thermopile embedded in a metal block behind 
the receiving surface. Similarly we may write 


G( 28) oT + ky Tx 7.) (8) 
From equations (7) and (8) it follows: 
G(1s) — G(2s) = — + — (9) 


For small temperature differences, the first right-hand term 
can be linearized and combined with the second term to give 
G(1s) — G(2s) = ko — T2,) (10) 


A similar equation holds for the condition that the reference 
surface is in the test position. Therefore 


C, (Ti. Tx) 

= C,) — Tw) 

The galvanometer reading A is, because of the small tempera- 
ture differences involved, directly proportional to the temperature 


(11) 
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difference between the hot junction (7), or 72,) and the cold junc- 
tion 7',, which is maintained at a constant temperature. There- 
fore 
A(is) = — T.) (12) 


Similar expressions may be written 
This leads to the final equation 


for the other readings. 


C, A(is) 


A(2s) 
C, A(ir) 


A — A(2r) 


(13) 
This equation may be simplified further since the reference 
material was selected to have similar directional reflection charac- 
teristics asthe sample. ThereforeC, = C,. The terms A(2s) and 
A(2r) in equation (13) have been found to be negligible since the 
test surface and thérmopile cold junction were kept at the same 
temperature by the cooling water. The reflectivity was then cal- 
culated from the equations 


(14) 


Surface Transmissivity. The data from the experimental proce- 
dure for transmitted energy can be analyzed in a manner similar 
to that for the reflected energy. Performing heat balances on 
the thermopile for both the incident i and transmitted r+ 
energies and relating these to the thermopile output A result in 
the following expression for the transmissivity of the test surface: 


A(r) 
7 


where A(r) is the galvanometer reading with transmitted energy 
impinging on the thermopile and A(z) is the galvanometer reading 
with the total incident energy impinging on the thermopile. 


DISCUSSION 


M. N. Aref* 


The authors are to be congratulated for contributing additional 
information on the radiation properties of porous fabrics. It 
appears from the comparison technique adopted by the authors 
that, in order to arrive at the absorptivity and hence the emis- 
sivity, direct measurements were obtained first on the reflectivity 
and transmissivity. The difference of both factors from unity de- 
termined the absorptivity. In this method, however, any errors 
involved in the calculation of the absorptivity will be due to the 
measurements of both reflectivity and transmissivity. It would 
be helpful to know the magnitude of such errors, which might be 
attributed to a possible heat loss through the radiometer aper- 
ture and an absorbed heat by the comparatively wide radiometer 
mirrored surface. 

As is well known, however, the emissivity of a surface may not 
be equal to its absorptivity, unless its monochromatic absorp- 
tivity is considered independent of both radiation wave length 
and temperature. On the other hand, a direct measurement of 
the emissivity of such surfaces could be determined by either 
one of the two methods presented in an earlier paper by the 
writer,® with a slight modification of the apparatus for the 
normal solar radiation measurement. This is basically ac- 
complished by placing a solar radiation or a black-body radiation 
source at the same temperature as the sample, at a fixed distance 
from the measuring radiometer, provided that both source and 
sample bear the same configuration to the radiometer. Both 

4 Research Department, Foster Wheeler Corporation, New York, 
N.Y. Assoc. Mem. ASME. 

5 M.N. Aref, “Emissivity Measurements of Industrial Surfaces Due 


to Thermal Radiation,” Paper No. 58—HT-18, presented at ASME- 
AIChE Heat Transfer Conference, Chicago, Ill., August, 1958. 
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radiations from source and sample could then be directed at the 
radiometer separately and their ratio would render the emissivity 
directly, thus eliminating the sources of error in the measure- 
ments of the reflectivity and transmissivity. 


Authors’ Closure 


The first two questions asked by Mr. Aref can be answered 
directly from the Analysis of Data section. The amount of 
energy lost through the aperture is taken into account by the fac- 
tor C (amount of energy reaching thermopile from test surface 
divided by actual amount leaving test surface). Referring to 
previous work of the authors, Ref. [3] of the paper, the value of 
energy lost through the aperture for various surfaces varies from 
1 to 4 per cent. The amount of energy absorbed by the hemi- 
spherical surface is also taken into account in the analysis by the 
factor (1 — p,,). 

The authors are well aware of available methods for directly 
measuring the total emissivity and indeed have made such 
measurements in the past.67 However, if a radiometer is sighted 
on @ porous surface such as a parachute fabric, the energy reach- 
ing the radiometer comes not only from the porous sample but 
additionally some energy is transmitted through the sample; it 
is extremely difficult to separate the emitted energy from the 
transmitted energy and, consequently, the authors chose to de- 
termine the absorptivity as outlined in the paper. Mr. Aref ap- 
pears to have missed the point that the resulting measured ab- 
sorptivity for infrared radiation coming from a black body at 350 
F does indeed yield the emissivity if the monochromatic absorp- 
tivity is invariant with temperature (it can still vary with wave 
length). The proof is straightforward: 

* E. R. G. Eckert, J. P. Hartnett, and T. F. Irvine, Jr., ‘‘Measure- 
ment of Total Emissivity of Porous Materials in Use for Transpira- 
tion Cooling,”’ Jet Propulsion, vol. 26, 1956, p. 280. 

7 T. F. Irvine, Jr., J. P. Hartnett, and E. R. G. Eckert, ‘Solar 
Collector Surfaces With Wavelength Selective Radiation Character- 
istics,” Solar Energy Journal, vol. I1, July-Oct., 1958, p. 12. 


Kirchhoff’s law states that the spectral emissivity of a surface at 
a temperature 7’ is equal to the spectral absorptivity at this same 
temperature 


= Ayr (15) 


The definition of the total absorptivity @ of a surface at a tem- 
perature 7’, for incoming radiation J): 


Qy 
f JindX 
0 


The definition of the total emissivity of a surface at temperature 


T is given by 
TON 


en 


f, eon, TAX 


€bx,7 = amount of energy emitted by a black body at tem- 
perature 7’ at wave length A, Btu/hr ft? 


a= (16) 


where 


Now in the infrared experiments described herein the incoming 
radiation J; corresponded to black body radiation at 350 deg 


Jin = br, 350 F (18) 


The surface temperature of the parachute material was ap- 
proximately 100 F and, if it is assumed that the spectral absorp- 
tivity is independent of temperature, 


@,100 F = @,350 F = €\,350 F (19) 


Substituting (18) and (19) into equations (16) and comparing 
with equation (17) it is readily seen that the measured absorp- 
tivity directly yields the emissivity. 
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An Experimental and Analytical Study of 
Vortex-Flow Temperature Separation by 


Superposition of Spiral and Axial Flows 
Part 1 


This paper reports on an experimental and analytical study of compressible fiow in a 
uniflow vortex tube. Part 1 deals with an experimental study, Part 2 with the analyti- 
cal study. Its purpose is to provide a better understanding of the separation of a gas 
stream into regions of high and low stagnation temperatures, there being at present 
little agreement as to the theory of operation. The problem is first approached from 
the experimental standpoint. <A large, multipurpose vortex tube is so designed and . 
built that pressure, temperature, and velocity traverses can be taken at six different 
stations throughout the length of the tube. Pressure, temperature, and velocity traverses 
are taken by means of hypodermic probes. Velocities are checked by means of a minia 
ture hot-wire anemometer. Data are taken for different runs of inlet pressures and ° 
plotted against radial distance. Flow visualization 1s obtained by means of liquid in- 
jection. The analytical study consists of using superposition for the solution of the 
flow equations. It begins with potential vorter flow in the plane. The solution of 
this flow is characterized by the existence of sonic or limit circles. Superposition of a 
sink flow to the vortex solution yields a spiral flow in the plane. The general solution in 
space is obtained by addition of a uniform axial velocity to the spiral flow. When 
viscosity effects are considered, the potential vortex changes into a forced vortex, and the 
solution becomes a superposition of a viscous compressible sink to a forced vortex. 
Performance or stagnation temperature separation is expressed as function of the ratio 
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Introduction 


HE PURPOSE of this work is to provide a better 
understanding of the separation of a gas stream into regions of 
high and low stagnation temperatures in a vortex field. Its 
motivation is the perplexing phenomenon known as the Ranque- 
Hilsch effect, which occurs when a gas expands in a centrifugal 
field such as exists in cyclone separators. This effect is a fairly 
recent discovery, having first been reported in 1931, and redis- 
covered in 1946. Its history follows the familiar pattern of initial 
As a 
phenomenon, relatively little is known concerning it, except for 
its spectacular effect of producing hot and cold air simultane- 


enthusiasm, then apathy, and later, renewed interest. 


' Paper is based on some of the experimental work done by the 


author at the University of Michigan. 

Contributed by the Heat Transfer Division and presented at the 
Semi-Annual Meeting, Detroit, Mich., June 15-19, 1958, of THe 
American Socrery oF MECHANICAL ENGINEERS. 

Nove 


Statements and opinions advanced in papers are to be 


understood as individual expressions of their authors and not those 
Manuscript received at ASME Headquarters, May 
Paper No, 58 


of the So« ety 


14, 1958 SA-71. 


Nomenclature 


of vortex strength to sink strength. 


ously. 


Despite various hypotheses advanced, there is to date, 
no general agreement as to its theory of operation [12].2 Conse- 
quently, there is a need for a systematic study of the phenomenon 
from the ground up. The present paper represents such a study. 
It is hoped that it will be of help in the understanding of the 
Ranque-Hilsch effect and of vorticity in general. 

The vortex tube is a very simple device which has no moving 
parts, but merely consists of a straight length of tubing with a 
tangential entry for the supply air as shown in Fig. 1. The enter- 
ing stream of compressed air separates into two streams of dif- 
ferent temperature; the hot stream leaving through the periph- 
ery, and the cold stream leaving through the center. Except 
for its spectacular effect of producing hot and cold air simul- 
taneously, relatively little is known concerning its operation. 
So far, analyses have been based on very little systematic experi- 
mentation because of the difficulties involved in gathering data. 
Also, the equations describing the flow are difficult to handle, for 
they are of the spatial viscous compressible-flow type. Never- 
theless, since 1946, when it was brought over from Germany, 


? Numbers in brackets designate Bibliography at end of Part 2. 


A = area ¢ = speed of sound Ir, ¥, 2 = co-ordinate axes 
AK = constant in irrotational flow m = mass Y = ratio of specific heats 
p = pressure = circulation 
: q = velocity in hodograph plane 6 = velocity angle in hodograph 
= strength of source or sink . 
r = radial distance plane 
= gus ¢ sts 
KR = gas constant t = time = density 
S = entropy u,v, w = components of velocity in physi- velocity potential 


absolute temperature cal plane 
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= streamline function 
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writers of widely different backgrounds have been interested in 
the phenomenon [5, 6, 8, 17, 22, 30. 31] 

The history of the vortex tube dates back to 1931, when it was 
first reported by G. Ranque in a paper read to the Société 
Frangaise de Physique. He was a metallurgist at a steel works 
in Montlugon, Central France, and he had noticed the vortex 
cooling effect in connection with cyclone separators. He con- 
structed a device to duplicate the effect, and applied for a French 
patent in December, 1931. A similar patent was applied for in 
the United States in December, 1932, and was granted in March, 
1934. Subsequent developments, however, brought the disap- 
pointing realization that the vortex tube was inefficient as a 
refrigerator, for nothing more was heard of it until 1946, when 
R. Hilsch, of the Physikalischen Institute, Erlangen, Germany, 
published a paper relating his experience with the device [16]. 

In June, 1945, at the end of the war in Europe, C. W. Hansell 
investigator for the United States and British Technieal Indus- 
trial Intelligence Committee, visited Erlangen University and re- 
At about the same time, R. M. 
Milton, of Johns Hopkins University also visited Erlangen, and 
brought back a model and a thesis by Hilsch. The thesis was 
translated by I. Estermann, of the Carnegie Institute of Teeh- 
nology, and the translation was circulated through Wright Field. 
Hilsch’s original paper gave some performance data and optimum 
dimensions for two vortex tubes. 


discovered the vortex tube. 


Widespread American interest 
has been given the vortex tube ever since Milton published a 
short descriptive article in May 
1946. The extreme simplicity of the device and the spectacular 
separation of hot and cold streams suggested that it might re- 


Refrigerating Engineering, 


How- 
ever, subsequent investigations showed that the power required 


place some of the more complicated refrigerator designs. 


to operate the vortex tube was many times that required by a 
conventional refrigerator. Nevertheless, other applications have 
been mentioned for the device, notably as regards jet-aircraft 
{22}. Thus, interest in the phe- 
nomenon continues, and in common with other efforts this paper 
aims to clarify the Ranque-Hilsch effect and to extend the knowl- 
edge of vorticity in general. 


Experimental Program 


Vortex-Tube Design. 


cooling and thermometry 


The main consideration in the test program 
was the design of a vortex tube of sufficient size and flexibility 
as to allow velocity, temperature, and pressure traverses to be 
taken without causing major disturbances in the flow field. The 
inherent limitations to the experimental study of vortex flow, 
until recently, have kept the systematic investigation of the vor- 
Both Ranque’s and Hilsch’s 
original models were of very small diameters (4 to 18-mm tubes 
with 2 to 7-mm orifices) wherein fairly impressive effects could 


tex tube from progressing very far. 


be obtained with relatively low or moderate supply pressures. 
Such small-size models, however, are not suitable to any basic 
study of the vortex phenomenon, since they do not lend them- 
selves to any velocity, pressure, or temperature traverses. To 
perform these operations, considerably larger size models have 
to be designed, and this in turn requires high flow rates of supply 
air. In addition, great care must be given to the fact that what- 
ever probes are inserted in the tube, the flow pattern must not 
be disturbed. 
With these considerations in mind, it was decided to design 
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Vortex or Ranque-Hilsch tube 


2-in-diam vortex tube made of lucite. This is considerably 
larger than the majority of previously reported designs. The 
choice of lucite, of course, is to permit flow-visualization studies 
in addition to the systematic taking of traverses. The installa- 
tion evolved is flexible for many kinds of measurements at dif- 
ferent stations along the length of the tube. It is shown in Figs. 
2,3, and 4. 

Compressed air at room temperature is fed tangentially into a 
center block, Fig. 3, from which it spirals downward along the 
vortex tube, Fig. 2. The vortex tube is so designed that a num- 
ber of sections are fitted together to form a tube of any desired 
length. Stagnation stagnation static 
pressure, and velocity traverses are taken at stations along the 
The exit end of the tube is fitted with a 
cone-shaped valve, Fig. 4, which is movable in and out to regulate 
the flow in the hot tube. 

Test Installation. Referring to Figs. 6 and 7, compressed air is 
supplied by means of valve no. 1 to a New Jersey constant-head 
It then passes through a DeVilbiss regulator and 
drier, and is fed by a flexible hose to the center block of the vortex 
tube. Pressure gage no. 1 and thermocouple no. | are installed 
between the Jersey flowmeter and the DeVilbiss drier, and a 
mercury U-tube is installed at the outlet of the drier so as to record 
the press're of the air at entrance to the vortex tube. 


temperature, pressure, 


tube shown in Fig. 5. 


flowmeter. 


A bypass 
circuit, consisting of valve no. 3, pressure gage no. 3, is included 
for inlet pressures which exceed the capacity of the DeVilbiss drier. 

No further eq iipment than the foregoing is needed for the study 
of vortex flow. except that a miniature pitot tube, Fig. 10, is used 
The flow in the hot tube 
is obtained very simply by taking the difference between the 
supply flow and the flow in the cold tube. 

Instrumentetion. For the experimental work described in this 
paper, special velocity, pressure, and temperature probes were 
The key feature is a probe assembly, Fig. 8, which can 
he inserted at anv station along the length of the vortex tube, 


to measure velocities in the cold tube. 


designed. 


and constructed in such a way that a hypodermic needle may be 
raised, lowered, or even completely revolved within the flow field. 

The hypodermic needle is raised or lowered by means of a 
Brown and Sharpe 605 depth gage and a slider as shown. The 
probe assembly is so constructed that the hypodermic needle not 
on'y may be moved radially into the flow, but also may be re- 
volved completely so as to be sensible to the direction of flow. 
Fig. 8 shows the probe that is inserted at stations 1, 2, 3, and 
so forth, whereas Fig. 3 shows the probe that is inserted in the 
entry block. 

The velocity, pressure, and temperature probes are shown in 
Figs. 9 and 10. Velocities may be obtained by measuring the 
stagnation and static pressures, respectively, but it is more con- 
venient to use a hot-wire anemometer, the theory, electric circuit, 
and contro! panel of which are standard knowledge. 

The hot-wire anemometer itself is a miniature version of the 
Referring to Fig. 
9, it consists of inserting two insulated copper leads through a 
stainless-steel hypodermic tubing (no. 18 gage), and soldering 


more bulky types made by commercial houses. 


two sewing needles to these leads. The needles are held apart 
and in rigid position by binding them with fine thread and coat- 
ing the whole with lavers of electrical insulation which are 
allowed to dry and harden. The hardened base (shown as a 


small white mass in Fig. 9) is reduced to proper size by means of a 
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Fig. 2 Vortex-tube design 


CENTER BLOCK 
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VORTEX TUBE 
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— 


STATIONS 


Fig. 4 Tail-cone assembly 


grinding wheel, and finally it is given several coverings of air- 
Across the sewing needles 
(which can be seen protruding from the base) is soldered a fine 
Wollaston wire ranging in size from 0.00015 to 0.00030 in. diam. 
To obtain high sensitivity to velocities, the center portion of the 
wire is etched in a solution of nitric acid. 

The pressure probes, Fig. 10, consist of a static-pressure hypo- 


craft dope to smooth the surface. 
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COUNTERFLOM VORTEX 
Fig. 5 Traverse stations 


dermic steel tubing (no. 18 gage) polished and left open at the 
end, and of a stagnation-pressure hypodermic steel tubing (also 
no. 18 gage), whose end is soldered closed and square cut, but 
having a small hole which may be directed into the line of flow. 
By revolving the stagnation tubing on its stand, the impact pres- 
sure yields a maximum reading, which indicates the direction of 
flow. 

As for the stagnation-temperature probe, it consists of a stain- 
less-steel hypodermic of the same size (no. 18 gage) as those used 
for velocity and pressure measurements. Two dissimilar but 
insulated leads (iron and constantan, no. 36 gage) are inserted 
through it. The ends of the leads are bared and fused together 
by acetylene torch of carbon arc, then pulled back close to the 
end of the tube to make the assembly compact. To keep dis- 
turbance of the flow field at a minimum, no further elaboration 
is included in the design of the probes. All probes, however, are 
calibrated before use, the velocity probes in a wind tunnel, and 
the temperature probes by comparison with a probe of unity 
recovery factor. 

Test Configuration. The primary objective of the experimental 
program was to obtain basic knowledge of the velocity, pressure, 
and temperature distributions in the vortex field of flow. The 
emphasis of the study was directed to a vortex which is generated 
at the entrance block, and which proceeds in one direction down 
the tube to be discharged through the exit at the far end of the 
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Fig. 7 Flow circuit 


tube. This is the so-called uniflow type of vortex tube which 
heretofore has not been reported upon extensively in the litera- 
ture. The design of the tube, however, is such that it can be 
converted easily into a counterflow vortex tube, as demonstrated 
in Figs. 1 and 2. This is the commonly known Ranque-Hilsch 
tube. 

To cover adequately most of the geometry of flow, the tests 
were subdivided into flow configurations based on the selected 
length. Each length is denoted by a letter. The configurations, 
labelled A, B, C, and so on, are shown in Fig. 11. 

Within each geometrical configuration, traverses were obtained 
for inlet pressures of 10, 20, 30, and so on, psig, maintained by the 
pressure-regulator valve. Of the measurements gathered, the 
results for the 10 and 20-psig runs are chosen as representative, 
and they are reproduced for the greater part in this paper. 

Test Procedure. = Hach run consisted of the following operational 
procedures: a) The valves of the flow circuit, Fig. 7, were 
cracked open and the inlet valve to the vortex tube adjusted for 
the desired pressure of the run. 6) The exit cone was turned all 
the way in, then backed out fairly slowly until the greatest 
Ranque-Hilsch effect is obtained and left in that position for the 
remainder of the run. c) After allowing time for steady-state 
conditions to be reached, the process of traversing the tube was 
undertaken. Beginning with the probe assembly, Fig. 8, placed 
in position at station 1, the hypodermic probes for velocity, pres- 
sure, and temperature were introduced, one at a time, into the 
vortex tube and positioned at every tenth of an inch along the 
radius. The reason for introducing only one probe at a time 
was to keep the disturbance of the flow field toa minimum. For 
each radial position of a probe, two readings were taken, one 
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corresponding to the probe on its way from the tube wall to the 
center ¢f the tube, and the other corresponding to the probe on 
its way back from the tube center to the tube wall. The two 
readings were averaged, and this is the value that was recorded 
as being the probe reading at a given radial position. Since the 
vortex tube had a radius of 1 in., 20 readings were thus taken per 
item, and 10 averages were recorded for that item in the table 
of results. With the existence of the three items of velocity, 
pressure, and temperature, this meant 60 readings per traverse 
per station. With 6 stations along the length of the vortex tube, 
a given run involved no less than 360 readings or a 180 averages 
to be entered in the tabulation of results. With the additional 
factors that time had to be allowed for readings to stabilize and 
ambient conditions to be fairly uniform throughout, it can be 
seen that it is not easy to make a good run. At that, decent runs, 
corresponding to inlet pressures of 10, 20, and 30 psig were ob- 
tained. The results of the 10-psig run have been plotted in the 
form of curves in Figs. 12 to 17, and the remainder left in tabu- 
lar form as Tables 1 and 2 for the sake of economy. 

No runs were made above 35 psig because of the large con- 
sumption of air and of the noise involved. It should be noted 
that the readings near the center of the tube and near the wall of 
the tube were very erratic, and for that reason, neither traverse 
curves nor traverse data are drawn or given for radial distances 
less than 0.3 in. or more than 0.9 in. 

Axial Velocity. ©The velocities up to this point represent the 
magnitude of the velocity vector. The velocity probe, however, 
can detect direction. It can be rotated 360 deg on its stand, 
Fig. 8, so that in its position of maximum reading, the direc- 
tion perpendicular to the hot wire indicates the direction of the 
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Fig. 9 Velocity probe 
velocity. With this information, it is possible to resolve velocities 
flow 


station; 


into axial and tangential components. Fig. 19 shows the 


angle and the axial-velocity traverses at a typical 


namely, 2 

Traverses at other stations show curves of similar shape, and 
hence are not drawn. For convenience, the stagger or tangential 
Inci- 
dentally, it was found that the flow angle is rather independent 


angle (90-a@) is plotted, rather than the flow angle itself 


Fig. 


19 also shows that the axial velocities are uniformly small when 


of changes in inlet pressure over the range of values tested 


compared to tangential velocities or to total velocities as cross 
plotted in Fig. 18 


velocities are 


As for distribution, the predominant axial 
concentrated in a small annular region near the 


wall of the vortex tube. The axial velocity falls off sharply 
toward the center of the tube, and actually indicates a reverse 
flow near the center. This explains partly how the hot stream 
travels in one direction, and the cold stream in the other in the 
case of previous investigations of the vortex tube 

A reasoning concerning the wall friction throws considerable 


jight on the tangential-velocity distribution. 
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Fig. 10 Temperature and pressure probes 
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data reveals that the tangential or total-velocity peaks do not 
coincide with the axial-velocity peaks, but are located more 
toward the center of the tube. This can be explained on the 
basis that wall friction causes a decrease in both the axial and 
tangential velocities near the wall as the flow moves along. 
From continuity considerations, this results in an increase in axial 
and tangential velocities in regions farther away from the wall. 
It indicates that some air must move in toward the center as 
In so doing, the air tends to 


conserve its angular momentum, thereby resulting in an increase 
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Fig. 18 Axial variation of pressure, temperature, velocity; configuration A, Pi... = 10 psig 


static pressure in the neighborhood of the center line must in- 
crease from station | to station 6. Such a pressure gradient ex- 
plains the tendency for backward flow near the center of the tube, 
especially in the case of the counterflow vortex tube [16, 25]. 

Nor is the foregoing the only case of data cross checking. The 
axial-velocity distribution may be used in conjunction with the 
local values of the density to determine by integration, the mass 
rate of flow down the tube. This mass rate can then be com- 
pared with the mass rate obtained from the flowmeter. Simi- 
larly, an energy balance may be performed by taking the product 
of the meter flow rate, the specifie heat, and the inlet tempera- 
ture, and comparing with the integrated energy flow at various 
cross sections along the tube length. This was done on several 
occasions, and the results agreed well enough within experimental 
limits. This is one of the many advantages occurring from the 
“{nternal”’ approach to the study of the vortex tube. 

Flow Visualization. Construction of a vortex tube of lucite 
afforded good opportunity for the study of flow patterns by 
means of visual methods. Accordingly, a smoke chamber was 
first constructed and connected to the inlet hose. The injected 
smoke showed up well at low speeds, but at moderate or high 
eabiee brated speeds it became invisible. Also, there was considerable dif- 
ficulty in getting the smoke into the vortex-tube center block. 
The pressure in the smoke chamber being low compared to the 
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Fig. 19 Flow angle; axial velocity 


in angular velocity in those regions farther away from the wall. 
The tangential velocity, however, drops off at the center of the 
tube, for as will be shown in Part 2, there is theoretically no flow 
“in the plane,”’ but only axial flow. 

Interrelationship of Internal Data. ‘The experimental study of the 
vortex tube “from the inside,”’ as is the case in the present in- 
vestigation, affords some interesting cross checking of data. It 
was pointed out previously that the plot of axial velocities 
showed some negative values in the neighborhood of ‘the axis, 
indicating the existence of a backward flow. This may be veri- 
fied by considering the pressure relation for circulatory flow: 


(1) 


Figs. 12 through 18 show an increase in tangential velocities near 
the center of the tube. Since Equation (1) must be satisfied, it 
follows that there must be an increase in the pressure differential 
between the center line and the region immediately away from it 
as flow progresses down the tube. Now, static-pressure meas- 
urements along the center line itself showed fairly constant 


values along the length of the vortex tube. Therefore, the Fig. 20 Flow pattern upon entry 
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Fig. 21 Limit circle-counterfilow vortex tube 


Fig. 22 Limit circle-water injection 
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air-supply pressure and to the vortex-tube pressure, the smoke 
had a tendency to flow out of the tube rather than into the tube. 

After smoke, balsa sawdust was tried but it did not show up 
well at high speeds; also, it covered everything in the vicinity 
with a coat of dust. Next, fine confetti was tried, with a laundry 
bag covering the vortex-tube exit cone to prevent the confetti 
from spreading to the adjacent laboratory space. It worked well 
for short bursts of speeds, but then the confetti, in rubbing 
against the lucite, developed static electricity, and it simply 
clung to the walls of the tube to no avail. 

The final and most successful method of flow visualization was 
The 
injections were made by means of «a horse syringe, and clear 
water was first kept the installation clean. 
Later on, for the sake of taking photographs, colored liquid and 
milk were used 


the injection of clear or colored water into the vortex tube. 
used, beeause it 


This entailed washing and cleaning the entire 
installation after an injection, however, it was not too high a 
price to pay for the acquisition of good photographs. Fig. 20 
shows the flow pattern immediately after entry into the vortex 
tube. The circle in the field of flow is the limit circle, and repre- 
sents the inner boundary of the vortex flow. Its existence is a 
characteristic of compressible flow, and its theory will be de- 
veloped in Part 2 of this paper. Fig. 21 shows the same condi- 
These 
are interesting pictures, for the existence of limit circles (circles 
beyond which the flow does not penetrate) has not 


tion, but for the case of the counterflow vortex tube. 


been dis- 


cussed in previous investigations. Fig. 22 shows another ex- 


ample of limit line in plane flow. However, it was obtained by 
injecting a trace of clear water, and it represents one of the few 
cases where a ‘“‘clear’’ trace photographed well. 

Fig. 25 shows a streamline in three-dimensional flow. It was 
obtained by means of clear liquid, but using a short configura- 
tion; namely, configuration D. A similar streamline in space, 
but using a longer configuration will be shown in Part 2, where it 
is included to correlate the theoretical study with some of the 
experimental results. 


Journal of Heat Transfer 


Conclusions 

While the interpretation of results is deferred to Part 2 of the 
paper, the conclusions gathered from the experimental study 
may be stated as follows: 


1 Although it is difficult to take measurements inside a vortex 
tube, traverses are possible if a big enough tube is designed. 

2 Hypodermic probes are ideally suited for the purpose of 
taking measurements inside the tube, but probes must be used 
separately so as not to disturb the flow. This, however, means 
numerous changes during a particular run and makes the experi- 
mental work tedious. 

3 It is hoped that future investigators will design a multi- 
purpose probe, which can be used to measure temperature, pres- 
sure, and velocity without interchanges, and yet be small enough 
not to alter the flow field. 

t Llow visualization is possible by means of liquid injection 
An envelope streamline may be obtained by liquid injection, but 
other streamlines are difficult to obtain either by liquid injection 
or smoke injection. Liquid injection shows that vortex-sink flow 
is formed in the entry block rather than farther down the tube 

5 Results of the experimental study are in general agreement 
with past experimental work that has been done in this field; 
Further 
work, both experimental and analytical, will add much to the 
understanding of vorticity temperature separation. 


however, no claim is made that the study is exhaustive. 


DISCUSSION 
Frank Kreith® 


The experimental phase of this investigation on the separation 
of a gas stream into regions of high and low stagnation tempera- 
tures adds valuable data on the macroscopic characteristics of a 
uniflow vortex tube. Those interested in phenomena related to 
energy separation in vortex-type flows are undoubtedly looking 
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forward to the presentation and publication of the theoretical 
phase of this work without which it is difficult to evaluate the 
over-all significance of the author's contribution, and to comment 
with confidence on the experimental data presented in this paper. 

The macroscopic mechanism of energy transfer has been quali- 
tatively explained by Schultz-Grunow [1],‘ Fulton [2], Van 
Deemter [3], and most recently by Deissler [4]. Since, except 
at very low Reynolds numbers, the energy separation in laminar 
flow appears to be quite small, interest is focused on turbu- 
Sub- 
stantially, two mechanisms seem to control the distribution of 
total energy in turbulent vortex flow. There is first of all tur- 
bulent mixing in the presence of a centrifugal pressure field 
which tends to build up an adiabatic temperature distribution. 
But near the inlet nozzle, as a result of friction forces acting be- 
tween faster moving gas in the center and slower moving layers 
of gas at the periphery, there occurs radially outward a transfer 
of energy as work. 


lent flow where the energy separation may be appreciable. 


The internal energy in the outer region in- 
creases therefore at the expense of the energy stored in the inner 
core, Simultaneously, however, the inner core is heated by con- 
In other words, the tem- 
perature profile is determined by the ratio of work flux to heat 
flux. 

The validity of this model can, of course, only be established 
by comparing predictions based upon it with experimental data. 
Unfortunately, in order to calculate the temperature distribution 
it is necessary to make additional assumptions and simplifica- 
tions. 


ductive and convective heat transfer. 


One must, first of all, assume a velocity profile which for 
turbulent vortex flow can only be based on very secant experi- 
These data are indeed a welcome addition 
But in addition to the velocity profile one must also postu- 
late some sort of model of the turbulent heat-transfer mechanism 
in the presence of a centrifugal force field. 


mental evidence 
here. 


As shown in reference 
[5] there are at least three transferable properties which can be 


* Numbers in brackets designate References at end of discussion. 
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used to evaluate a coefficient of eddy conductivity. The most use- 
ful of them appears to be the so-called forced vortex parameter. 
But even if one accepts the forced vortex parameter and neglects 
the peculiarities of boundary layers in swirling flow [6], there 
still remain a number of puzzling effects, such as the qualitative 
difference in behavior between compressible and incompressible 
flows in a diabatie vortex flow [7], which can only be partially 
explained in terms of existing information on the energy transfer 
in vortex flow. It seems therefore important to supplement the 
author’s macroscopic measurements of time average values of 
pressure, temperature, and velocity distributions, by measure- 
ments of the fluctuations of these quantities. In particular, a 
knowledge of the effect of swirling on the radial velocity fluctua- 
tions would be imperative for a more complete analysis of the 
The approach used by Yeh [6] 
and Einstein and Li [8] may be applicable and should be explored 
further. 


energy-transfer mechanism. 
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Part 2' 


Part 2 reports on the analytical study. 


Vortex-Flow Temperature Separation by 
Superposition of Spiral and Axial Flow 


The free vortex motion of the gas upon 


entrance to the tube is mathematically superposed to a compressible sink to give a spiral 


flow in the plane. 


The characteristic existence of limit circles is corroborated by the 
experimental flow visualization. 
uniform axial velocity to the spiral flow. 
vortex is shown to change into a forced vortex. 
kinetic energy and maximum entropy. 


The solution in space is obtained by addition of a 
When viscosity effects are considered, the free 
The latter flow is one of minimum 
Energy considerations enable the determina- 


tion of an optimum cold air radius to give largest stagnation temperature separation. 
Significantly, this was the radius that gave best performance in the experimental 


program. 


Introduction 


Parr 1 reported on the experimental study of the 
vortex effect, and presented the test results. Part 2, which is 
the present paper, reports on the analytical study, and presents 
a mathematical solution for the differential equations governing 
the flow in the vortex tube. For greater understanding, it is 
recommended that Part 2 be kept in mind when reading the pres- 
ent paper. 


Analytical Study 


Consider the flow of air upon its entrance into the vortex tube. 


Its motion is as shown in Fig. 1. The analytical treatment of 


the problem of solving the flow equations will proceed along the 
following line of thought: The intrinsic equations of compressible- 
fluid flow, if tackled head-on in their general form, namely, three- 
dimensional, with viscosity effects included, involve mathemati- 
cal difficulties which prevent a direct solution. It is, therefore, 
only natural that simpler solutions be sought from simpler models 


' Part 1 was presented at the ASME Semi-Annual Meeting, De- 
troit, Mich., June 15-19, 1958. It reported on the experimental 
study as based on work done at the University of Michigan. Part 
2 reports on the analytical study. 

Contributed by the Heat Transfer Division and presented at the 
Annual Meeting, New York, N. Y., November 30-December 5, 1958, 
of THe AMERICAN Society OF MECHANICAL ENGINEERS. 

Note: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those 
of the Society. Manuscript received at ASME Headquarters, August 
7, 1958. Paper No. 58—-A-90. 


————— Nomenclature 


of flow, and by superposing the individual solutions, to arrive at 
the general solution. Using Prandtl’s concept of boundary layer, 
it is possible to initially ignore friction in the flow field outside the 
boundary layer, and to achieve exact solution based on potential 
flow. 

Free Vortex. Consider a perfect fluid moving in a circular 
path with no external torque applied: 


=m 
dt 


A / ig 
fRIBUTION 
17 


ot! 


(a) 


Fig. 1 Flow in circular path 


= area 
constant in irrotational flow mass 
Mach number = pressure 
strength of source or sink 
= gas constant tine 
= entropy 


absolute temperature 


Journal of Heat Transfer 


= speed of sound 
velocity in hodograph plane 
radial distance 


components 
physical plane 


co-ordinate axes 

ratio of specific heats 

circulation 

velocity angle in hodograph 
plane 

density 

velocity potential 

streamline function 


of velocity in 
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Fig. 2 Velocity, pressure, temperature distribution 
for free vortex 


Integration of Equation 


1) gives 
rv = K; v= (2) 


Equation (2) is characteristic of a free or potential vortex. The 


velocity distribution is of the hyperbolic type shown in Fig. 2(a). 
The pressure distribution is determined by considering the force 
balance in Fig. 1(4). The centrifugal force on the fluid element is 
balanced by the resultant force due to pressure; neglecting in- 
finitesimals of higher order than the first, the force balance in the 
radial direction is 

2 


y 
(dm) = 


(pdrdA) 


(dp (dA) = 


dp = 


Substitution of Equation (2 


dp = pK? ~ 4) 
ri 
For a perfect gas, p = RT’ and Equation (4) becomes 
‘ 


To integrate Equation (5), use is made of the fact that the total 


energy in a stream tube is equal to that in any other stream tube: 


(6) 


= 


dp K? dr 2e,K 


which is integrable by use of partial fractions, and yields a pres- 


sure distribution of the form shown in Fig. 2(4). As for the tem- 


perature distribution, it is obtained very simply from Equation 


(6), which vields 
2c, 


and is shown in Fig. 2(e). 
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Fig. 3 Circulation for free vortex 


“« 


Fig. 4 Fluid rotation at a point 


Circulation for Free Vortex. The potential, or free vortex just 
described, is a two-dimensional motion in which the streamlines 
are concentric circles, with the velocity along any streamline 
being inversely proportional to the radius. Referring to Fig. 3, 
the circulation, which is the sum of the products of elements di of 
a curve drawn in the plane of flow, and the corresponding com- 
ponent of velocity tangent to the curve, is 


r= (9) 


where p is the velocity. For a closed curve such as A-B-C-D, and 
noting that the line integral is zero along the sides AB and CD, 


the circulation is 


= (v + dv)(r + dr)dO — = divr\dé = 0 (10) 
since vr = K from Equation (2). The circulation for an element 
of area not enclosing the origin is thus zero. This result can be 
generalized to a finite area, since the line integral around the 
boundary curve of the finite area is the algebraic sum of all the line 
integrals around the elements of area comprising the finite area. 
The circulation is related to the fluid rotation in the manner 
shown in Fig. 4. Let the fluid rotation at a point be defined as 
the mean angular velocity of two mutually perpendicular fluid 
lines OA and OB drawn through P. During an infinitesimal 
time dt, OA rotates to OA’, and the vertical displacement is 


AA’ = dxdt 


or 


The angle AOA’, positive when measured counterclockwise, is 
then 


and its time rate of change or angular velocity is thus - Simi- 
Or 


-, so that the 


larly, the angular velocity of the fluid line OB is 


average rotation of the fluid element is 


ov ou 
= — 11) 
or oY 


Now, the circulation, as defined by Equation (4), 
pressed in terms of the velocity components u, r, «= 


rT = SF (udz + vdy) 


and this may be transformed by means of Green's Theorem to a 


surface integral as 
Ov ol 
- — —}drdy (13) 
A Or oy 


r = (udx vdy) = 
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into equation (3) gives 
dp K? dr 
= (5) 
RT 
Ov 
T+ = ¢, T + 
2 2r? 
“ 
Kk? 
2c,7 
Replacing 7’ by this value, Equation (5) becomes eC Sie 
2c r* 
7 
{ 
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A Fig. 5 Motion of fluid particle 
upon entrance to vortex tube 


Since the circulation for the free vortex, as given by Equation 
10), is zero, Equation (13) yields 


) 
= ( drdy; A = arbitrary 
or oY 


(14) 


Comparing Equations (11) and (14), it is seen that the fluid rota- 
tion at any arbitrary point upon entrance into the vortex tube is 
zero. In other words, the motion of fluid particles upon entrance 
to the vortex tube is as shown in Fig. 5, and the flow is potential. 
An exact solution for the flow pattern may therefore be expected. 

Potential Equation. The problem of determining the flow pat- 
tern for the vortex tube is solved by developing a differential 
equation of flow in terms of the velocity potential or the stream 
function which turns out to be linear. With the assumption of 
irrotational motion, the circulation, as calculated by Equation 
(10), is zero, and from this it follows that the line integral of the 
velocity between any two points of the flow field depends solely 
on the location of the points. Thus, o - di is the exact differential 
of a quantity @ such that in Cartesian co-ordinates 


(16) 
v-dl =d@ (15) 
(17) 
where u and v are the z and y components of the velocity. The 
equation of continuity, ; 


O(pu) O(pu) 
Or 


= 0 


may then be written with the aid of Equations (16) and (17) as 


fe) 
ag (b=) + 5 = + + + = 0 (18) 


Euler's equation, 


Op u? 
P dx + P = d ( 
Ox Oy 


becomes, in terms of the velocity potential, 


2 2 24 2 
dp = ) = ( (20) 


dp 


and this, with the expression for the velocity of sound c? ; 
dp 


may be rewritten as 
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from which 


= + (21) 
c 


p 
+ $,%,,) (22) 


Substitution of Equations (21) and (22) into Equation (18) 
gives 


(6,4, + $,%,) = 0 


or, upon simplification, 


(: ) + (: ) bw — = 0 (23) 


Equation (23) is the potential equation in terms of the velocity 
potential @; however, it is not linear, and will have to be trans- 
formed later into an equation which is linear. 

Stream Function Equation. Just as the condition of free vortex 
flow led to the existence of a potential function @, so can the 
equation of continuity be satisfied by the introduction of a stream 
function Y. This is because the equation of continuity 


(pu) + ( 0 
u 
Ox oY 


is satisfied if 


pu 


where y(z, y) is a point function. For, the order of differentia- 
tion being immaterial, 
oroy 
and it is obvious that Equation (26) satisfies the continuity equa- 
tion. Equations (24) and (25) define the stream function, but it 
is more convenient to alter them into 


where po is the stagnation density. The magnitude of the ve- 
locity at any point then becomes 


Just as a differential equation of flow was developed in terms 
of the velocity potential, so will a differential equation be de- 
veloped in terms of the stream function. First, the equation of 
continuity is automatically satisfied with the introduction of the 
function y, 
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A 
p + o,? 
dp = d\{ - 
c 2 
aed 
Pz 
- 
Ps Ov Ou 
ox oy 
oy 
u=- (24) 
oy a 
Or 
bs 
Pony 
u (<4) 
eo -F y, (28) 
p ox p 
+ y? 
p 
; 
; 


(29) 


oy 


(pu) + (pv) = 
or oy oyor 


Second, the condition of irrotationality 


oy or oy \p oY or 


becomes, after differentiation and rearrangement, 


or oy 


(30) 


Third, Euler's equation, 


1 2 2 
p 2 


dp . 
and the expression for the velocity of sound, c? = = yield, upon 


elimination of dp: 


dp = 


d(u? + = - (2) + v.'] 
2e 2c p p 


— + ¥,?) (31) 


Equation (31) enables the formation of the partial derivatives 
p, and p,, which upon substitution into Equation (30) results in 


p 


Equation (32), along with the relations between local and stagna- 
tion quantities 


(33) 


p 2 p 


completely determines the flow field. 
of a transformation. 

Hodograph Transformation. The exact solution applies to that 
region of flow outside the boundary layer. By excluding the lat- 
ter from the flow field, the behavior of the fluid becomes irrota- 
tional, and it is then possible to use superposition to obtain the 
general solution. Having obtained Equations (23) and (32), 
what technique may be employed for their solution? Mathe- 
matically speaking, Equations (23) and (32) are nonlinear, and 
there is consequently no direct solution. The situation, however, 
is not as hopeless as it seems. The equations are in terms of only 
two space variables, and for this condition, they are amenable to 
exact solution by means of the following transformation. 

Let z and y be the co-ordinates of the physical or actual flow 
plane, and q and @ be the polar co-ordinates of a velocity or 
hodograph plane. Then 


It will be solved by means 


oo (35) 


“= q cos = 
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(36) 


v =qsin@ = 
oy 


and the defining equations for velocity potential and stream func- 
tion: 


(37) 


give, upon differentiation with respect to z and y, 


= udr + vdy = q(cos Odn + sin Ody) 


dy = — p vdzr + q(—sin 6dz + cos Ody) 
Po 


p 
-udy = — 
Po 


whereupon 


q 

dy = + (39) 
q P 


Now, considering @ and W as functions of g and 6, the following 
may be written 


dod = dq + dé (40) 
dy = dq + (41) 


Substitution of Equations (40) and (41) into (38) and (39) gives 


oo ) 
dr = — —— de + dé 
q oq at 
Po dq + 40 
4q oq 


which, upon simplification, becomes 


po sin 6 
v,) dq 
q 4 


“os 6 sin 
4 < Po ys) d@ (42) 
q 
Similarly, there is obtained 
sin 6 cos 6 
dy = ( "“& + 2 v,) dq 
q p 4 
(= 6 po cos y ) (43) 
q p 4 
Equations (42) and (43) give 
cos po sin 0 cos 6 po sin 6 
= — r= — ve 
q | q p 4 
sin 0 py cos 6 sin 6 po cos 6 
Yq = Ye = + — 
qd q @ 


The objective is to obtain @ and y in terms of gand 6. To this 
end, z and y are eliminated in the foregoing equation by making 
use of the independence of order in partial differentiation. The 
result is 
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p or ius 
p 
4 |_| 
Po 
\ Pp 
A ; 
4 
Voy 
2 
y¥—1f po 
p dr = 
| 


po cos 0 
De 
p q 


q q 


po sin 


po cos 6 Yo (45) 
q Ldqg\p 


Multiplying Equation (44) by sin 0, Equation (45) by cos 6, and 


subtracting gives 
d 1 po 


Similarly, multiplying Equation (45) by cos 6, Equation (46) by 
sin 8, and adding gives 


sin 6 


de 
(46) 


(47) 


Equations (46) and (47) are the compressible-flow equivalents of 
the classical Cauchy-Rieman equations of incompressible flow. 


Next, using Euler's equation: : + qdq = 0, and the expression 


for the speed of sound: ¢? = , there is obtained 
dp 


d (= po dp po dp dp 
dq p 7 p dq - p dp dq 
po ( 1 ) p 
} (-—pq) = 


which, substituted into Equation (46) gives 


Po 


(1 — (49) 


q 
Finally, eliminating @ from this latter equation and Equation 


(47) by setting @,6 = o,, and making use of Equation (48) 
gives 


Po 


oe, = — (1 — M?)Qeo 
p 


q 
Po Po q 
or, upon simplification, 


+ + + (1 — M*) = 0 (50) 
Equation (50) is the hodograph equation for the flow in the vortex 
tube. It is a linear equation in the dependent variable w. 

Superposition. The general the solution of 
Equation (50) is as follows: (1) Assuming that a number of 
simple functions y,, Y2 have been found to satisfy Equation (50), 
more complex solutions are found by linear superposition: Y = 
cWi + (2) The partial derivatives and are obtained 
by partial differentiation of the solution Y = W(q, 8). (3) Using 
the equivalent Cauchy-Rieman equations (46) and (47), @, and 
ge are found. (4) Knowing @, and @¢, the partial derivatives 
ty, Te, and y,, ye of the physical variables are found from Equa- 
tions (42) and (43), and integrated to give the streamline. 


procedure for 
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Solution for Vortex Tube. The solution for the streamline 
pattern in the vortex tube is obtained by combining a vortex 
flow with a sink flow. The respective solutions of these in- 
dividual flows are well known, and are as follows: 

Free Vortex. A particular solution of Equation (50) is sought, 


such that depends only on g. If this is the case, then Ye = 0, 
and Equation (50) becomes 


which is integrated into 
y = CJ + Ce 
1 
partial derivative of x and y are 
sin 6 
q 
cos 6 sin 6 
w= 


q* q 


where J = 


cos 6 


| 


Integration of these gives the equations for the streamline in the 
physical plane: 


sin 6 


= ¢; 
q 


cos 6 


C2 
q 


Equation (53) shows that the streamlines of constant velocity 
are concentric circles in the physical plane. Choosing the con- 
stant to bec, = rq = r*c* where r* is the radius corresponding 


to a Mach number of unity, Equation (53) may be written as 


I sin 6 cos 6 
The velocity distribution, shown in Fig. 6, increases with de- 
dreasing radius to a maximum of g = 2.45 c*. This corresponds 
to a minimum radius, a zero density, and a zero temperature. 
In terms of vortex strength, the minimum or limit radius is 


The existence of a limit circle has been corroborated in Part 1. 

Sink. A solution for Equation (50) is sought, for which wy 
depends only on 8. Under these conditions, Y, = 0, and Equa- 
tion (50) becomes 


(55) 


the solution of which is 


(56) 


Fig. 6 Velocity distribution upon entrance to vortex tube 
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The streamlines of constant W are straight lines passing through 
the origin. The partial derivatives are 


po COs 4 ( 7h ) Po sin 6 
sin 
= Po ve = Po COs 
p p q 


They are integrated into 


po cos 

z= + ¢; 
@q 
po sin 

= ci 4 


The constants ¢, ¢;, and ¢, merely involve displacements of the 
entire streamline pattern, and without loss of generality can be 
set equal to zero. The parametric equations for the streamlines 
are thus 


Po cos 6 
| 


z= ¢, 


sin 0 


Equations (57) show, by division of y by x, that the lines of con- 
stant @ are radial lines through the origin. The velocity dis- 
tribution is obtained by choosing c; to satisfy the continuity 
equation, and plotting q/c* versus r/r* as shown in Fig.7. There 
are two branches to the curve, but only the subsonic branch is of 
interest, since the vortex tube was run at subsonic velocities. 
There is no value of g/c* forr < r*. In other words, the sonic 
radius is the limiting radius for the flow. 

Spiral Flow. The vortex flow and the sink flow are now com- 
bined to give a spiral flow in the plane. The result is 


Q 
¥ + J; R= 7 (58 
cos 6 r sin 6 
(p/po)(q/Co) 

(59) 

o sin cos@ | 


where T’ and Q are the strengths of the vortex and sink flow, re- 
spectively. A plot of the streamlines may be obtained by solving 
corresponding pairs of values of q/co and O for a given y, and then 
determining the values of z and y. The flow pattern in the plane 
is shown in Fig. 8. A similar picture was obtained in Part 1 by 
means of oil injection which left a visible trace on a lucite plate. 

Solution in Space. In the preceding sections, the equation 
of continuity, the equation of motion, and the relations for the 


existence of a velocity potential were combined to yield a dif- 


Fig. 7 Velocity distribution for sink 
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ferential equation of flow, which being in terms of two space 
variables, was solved by means of a hodograph transformation. 
For three-dimensional flow, similar technique does not bring 
success if the third space variable is brought in indiscrimi- 
nately. The procedure followed here [23]? consists in adding to 
the exact solution just derived an axial velocity such that it does 
not alter the flow equations in the plane. This velocity is of con- 
stant magnitude along the axis of the vortex tube and leaves the 
pressure, density, and temperature of the fluid unchanged. 

That the three-dimensional flow thus obtained obeys the same 
set of equations heretofore solved may be seen as follows: Let 
w = constant (in space and time) be the velocity component 
which is added to the two-dimensional flow. The question may 
be asked: Does the resultant motion give a flow pattern which 
is consistent with the flow pattern that has been solved in the 
plane? The answer lies with the equations governing steady- 
state compressible flow, namely, the equation of continuity: 


O( pv O(pw) 
+ = = (60) 
Or oY Oz 


Euler’s equations: 


Ou ou ou 1 Op 


u +w + — - (61) 


ow ow ow 1 cp 
u +w— + = () 635 
or oy oz p o 


and the equation of property change: 


Po 


For the case of two-dimensional flow, the last term in Equation 
(60) does not exist; for the superposed flow, it becomes 
O(pw) ow 


= - =0 


Oz Oz 


since p is independent of z, and w is constant along z. The con- 
tinuity equation is thus the same for both the two-dimensional 
flow and the superposed flow. 

For the case of two-dimensional flow, the third relation in 


Ou Ol 


Equation (63) does not exist, while the terms > in the first 


two relations are zero. For the superposed flow, the third rela- 
tion in Equation (63) vanishes due to the constancy of w with 


? Numbers in brackets designate Bibliography at end of paper. 


Fig. 8 Spiral flow in the plane 
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respect to z, y, z, and the independence of p with respect to z. 
Thus the two-dimensional flow and the superposed flow obey the 
same dynamical equations. 

Lastly, the property-change relation, Equation (64), remains 
unchanged, since p, p, and 7’ are not affected by the addition of a 
constant axial velocity, and since physical laws are not altered 
by uniformly translating co-ordinates. Consequently, the same 
equations describe both flows, and since these equations have 
been linearized and solved, what remains then is simply to super- 
pose the stream function of the two-dimensional flow with the 
stream function of the axial flow to obtain the general solution in 
space. This addition of a constant velocity in the axial direction 
does not alter the conditions of irrotationality. For, if the two- 
dimensional flow possesses a velocity potential @, then the addi- 
tion of a constant velocity w still leads to an irrotational flow field 
with @ + wz as the velocity potential. The plot of a streamline 
for the superposed flow is shown in Fig. 9. It is typical of the 
flow pattern that prevailed in Part 1 upon entrance into the 
vortex tube. 

Effect of Viscosity. The exact solution just developed is for 
that region of the tube (relatively near the entrance) where vis- 
cosity effects are negligible. The experimental results in Part 1 
show that, as flow proceeds down the tube, viscosity effects begin 
to take over, and the free vortex changes into a forced vortex. 
Using Kassner's concept of shear stress in circular turbulent flow 
[18], the following equation is written for the stress: - 


( dq q ) 
dr r 


where € is the turbulent exchange rate 


(65) 


Since the velocity dis- 


tribution upon entrance into the tube is q = 3 differentiation 
gives 
rdq + qdr = (): (66) 
and Equation (65) becomes 
(67) 
The shearing force and the moment of shearing force acting on an 


annular element of fluid a distance r from the center are, respec- 
tively ’ 


(68) 


; ( dq q 
oar r 
dr r 


The latter expression can be rewritten with the aid of Equation 
(66) a 


—4reK (69) 


2K 
M, = (- *) = 


Equation (69) indicates that an elemental moment of shear stress 
is constant and independent of r. 

Since there is no external torque applied to the flow upon en- 
trance to the vortex tube, the summation of internal elemental 
moments must be zero. Thus, 


f = 0 
0 


But M, is constant according to Equation (69); consequently, 
Equation (70) cannot be satisfied except for M; = 0, since no 
summation of quantities of like sign can be zero unless all the 
quantities are zero individually. In other words, the flow, while 
irrotational immediately after entry into the vortex tube, does 
not stay irrotational, but changes into something else so as to 
satisfy Equation (70), which is the (initial) governing relation. 
That this something else is rotational flow can be seen by applying 
the equations of rotational flow to Equation (70), and seeing that 
the latter is satisfied: Since the characteristic equation for rota- 
tional flow is g = wr, it follows that 


(70) 


(71) 


Substitution of Equation (71) into Equation (68) gives M,; = 0; 
consequently, Equation (70) is satisfied. Thus rotational flow is 
the stabler form of flow to which the initial flow transforms itself, 
even though both flows satisfy the conservation of energy princi- 
ple. It will be shown, however, that the rotational flow has the 
lesser amount of useful kinetic energy. In other words, its 
entropy is greater than that of the irrotational flow. The change 
from one flow to the other is due to viscosity, though it seems 
surprising that the relatively low viscosity of the air may have 
such an effect. The experimental results presented in Part 1, 
however, tend to verify this change, for the velocity traverses 
show a definite rotational distribution. 

Energy Distribution. The kinetic energy of each flow will now be 
computed. To this end, the velocities before and after conversion 
from a free vortex to a forced vortex are first obtained by equat- 
ing the moment of momentum of both flows to each other: 


K 
Minot = (r?)(2erpdr) = wpK(r2? — r,?) 
4 r 


1 


M vor -{ (w)(r?)(2mrpdr) = 


Setting Equations (72) and (73) equal to each other gives 


mpK(r.? — 


- Ge 


Fig. 9 Flow pattern in space 


aucust 1959 / 29 


: 
; 
“414 
.- 
ir r 
dr r — ri 4) (73) 
Da, = 
ng 
Journal of Basic Engineering 


whence 


2K 
o= (74) 
+ r;? 


Equation (74) is the relation between the strength of the forced 
vortex to that of the free vortex. The kinetic energy for the two 
modes of flow can now be obtained: 


= > (2mprdr) = In (75) 
© * 
™pw?* 
= 2larprdr) = (r2* — 
The latter may, by means of Equation (74), be rewritten as 
E, = (40) 


ro? + ri? 


Comparison of Equations (75) and (76) shows that the kinetic 
energy initially possessed by the free vortex is greater than that 
possessed by the rotational flow. The total energy of the two 
modes of flow, however, is the same, since there is no external 
transfer of heat or work from or to the system in both cases. 
Upon entry into the vortex tube, the fluid particles, because of the 
absence of external torque, move in concentric streamlines, with 
each streamline having the same amount of total energy. As the 
flow proceeds down the tube, viscosity effects begin to influence 
the flow (mainly because of the velocity increases to sonic speeds 
with decreasing radius), and the flow can no longer keep this con- 
stancy of energy distribution. The inner layers, having too much 
kinetic energy, must surrender part of their energy to the outer 
layers by being slowed down. The process of conversion ends 
with rotational flow being attained near the tube exit. The dif- 
ference in kinetic energy between the free vortex and the forced 
vortex is the energy taken up by friction to make the conversion 
from irrotational to rotational flow. 

The analysis, as given in the foregoing, is primarily mathemati- 
cal. It has the advantage of being easily derivable, and requiring 
no elaborate model for comparison. It is also corroborated by 
Eckert and 
rotational flow of the forced 
vortex type, both laminar and turbulent, produced tempera- 
ture separations that are independent of Prandtl numbers. 

Proceeding with the solution of the vortex tube, the problem of 
determining the flow pattern near the tube exit is obtained by 
adding an axial velocity to the rotational flow just as it was done 
for the irrotational flow. 


some case studies of boundary layers made by 
Hartnett [9], who showed that 


The plot of a streamline is similar to 
that for irrotational flow, except that, strictly speaking, a scalar 
velocity potential no longer exists, but each streamline is identified 
instead by a pseudo yp. 
Performance Calculation. 
not permit 


The present state of the field does 
a complete calculation of performance. The de- 
pendence on such factors as diameter of tube, length of tube, size 
of cold air orifice, exit throttling, angle of entry, ete., has not yet 
been systematically studied. It is therefore too early for the de- 
velopment of a single equation which would predict performance. 
The optimum size for the cold air orifice, however, can be de- 
termined from the present analysis. If Equations (72) and (73), 
which represent the moment of momentum for the irrotational 
and rotational flows, are equated, K, the constant of irrota- 
tional flow, replaced by rg, and w, the constant of rotational 
flow, replaced by q/r, there results 


WPT irror(T2? — ry?) = 


Grot(T2* — ri‘) 
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Fig. 10 Optimum size for cold orifice 


whence 


-T* 
= Girrot 
+ 


Equation (77) shows that as r 0, 0, and as r — r2, > 
24 irrot- 
irrotational to rotational, the velocity of the converted flow at 
the wall tends to twice the velocity of the otherwise irrotational 
flow at the same location. 


This means that as the flow changes in character from 


This, of course, is at the expense of 
the inner layers, whose velocities decrease from near sonic to close 
to zero. If 7; is the temperature of the air upon entrance to the 
vortex tube, then for the irrotational flow, the quantity: ¢,7; + 
K?*/2r,?, being the stagnation enthalpy near the wall, has the same 
value for all streamlines, and is shown by the horizontal line A-B 
in Fig. 10. For rotational flow, the quantity: ¢,7, + 2K?/r:*, 
while still the stagnation enthalpy near the wall, no longer is con- 
stant, but decreases from streamline to streamline as r tends to 
zero. Now, recall from earlier diggussion that the energy of 
rotational flow is equal to that of the irrotational flow minus the 
energy taken up by the viscous forces. Let the energy of irrota- 
tional flow, which is known from knowledge of the initial stagna- 
A — B — rs. 
Since the friction energy may be determined by taking the dif- 
ference between Equation (75) and Equation (76), the energy of 
rotational flow is known, and may be plotted as the area r, — C 
—- D The intersection of the two energy profiles gives a 
corresponding radius which is the optimum radius for largest 
temperature separation. This is the radius that was adopted for 
the experimental program, and significantly, it gave the best per- 
formance with respect to that variable. 


tion enthalpy, be represented by the area r; — 


The conclusions gathered from the present study on the vortex 
tube may be summarized by the following statements: 


(i) Upon entrance into the tube, the motion of the gas is that of 
a free vortex, because of the absence of external torque. 

(ii) The problem of the vortex tube is, to a great extent, 
amenable to mathematical solution through the superposition of 
a vortex and sink, followed by the addition of an axial velocity. 

(iii) The presence of viscosity changes the flow from a free to a 
forced vortex. The latter flow is one of minimum kinetie energy 
and maximum entropy. 

(iv) Though performance calculation for the general case is not 
possible at this stage, the optimum size for cold orifice may be ob- 
tained from energy considerations. 

(v) The report is a first study based on relatively detailed and 
accurate experimental data. In that respect, it has the ad- 
vantage of being free from conjectures in the analysis of the vortex 
phenomenon. 
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This paper and its earlier companion describing an experimental 
study até valuable contributions to an area of fluid dynamics that 
has been noted for its complexity and general muddiness. The 
treatment is broad and generally lucid, and it would appear that 
most of the effects having an important influence on the some- 
times peculiar characteristics of vortex flow have been taken into 
account. 

Two minor errors that do not appear in later equations and that 
do not affect any conclusions might be noted: v?/r? in Equation 
(3) should be v?/r, and K*/r in the equation following Equa- 
tion (7) should be K2/R. 

Will the author elaborate on details concerning the axial ve- 
locity distribution used? In a recent similar study‘ of vortex-flow 
energy separation, axial velocity was assumed constant with 
radius in a core region and, though of different magnitude, also 
constant with radius in the surrounding annular region. It also 
was shown that axial velocity varies linearly with length when 
tangential velocity is independent of length. In this paper, what 
is the form of the radial and axial variation of axial velocity? 

It also might be worth while for the author to outline the pres- 
ently unavailable knowledge which would be necessary for this 
theory to be numerically applicable to the design of vortex tubes, 
the prediction of their performance, and the correlation of vortex- 
tube data. 


F. Kreith’ 


The theoretical part of the author’s study of vortex-flow tem- 
perature separation deals with a problem which has been the sub- 
ject of a great deal of research in recent years. His paper pre- 
sents a welcome addition to the formal solutions available for 
calculating the potential flow pattern in a frictionless vortex flow. 
Unfortunately, the assumptions made in the course of the 
mathematical development are so restrictive that the final results 
are likely to deviate appreciably from actual conditions in a vor- 
tex tube, especially in turbulent flow which appears to be neces- 
sary for any substantial separation of hot and cold gas to occur. 

In order to ascertain the validity of the model postulated by the 


3 Long Range Planning Group, Oak Ridge National Laboratory, a 
Division of Union Carbide Corporation, Oak Ridge, Tenn. 

‘hk. G. Diessler and M. Perlmutter, ‘‘An Analysis of the Energy 
Separation in Laminar and Turbulent Compressible Vortex Flows," 
Heat Transfer and Fluid Mechanics Institute Preprints, 1958, pp. 
40-53. 

5 University of Colorado, Boulder, Colo. 
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author in his analysis, it would be of interest to know how closely 
the flow pattern observed visually in the vortex tube used in the 
experimental work (i.e., Part 1 of this study) agrees with predic- 
tions based on Equation (59) in this paper. In particular, since 
Equation (70) shows that the radial-velocity distribution given by 
the free-vortex hypothesis in Equations (1) and (2) cannot be main- 
tained in a real fluid, it would be important to know how it varies 
with axial distance from the entrance. It would also be of interest 
to know whether or not the experimental results show how to 
evaluate the strengths of the vortex flow and the sink flow in 
Equation (59). This information would enhance the practical 
value of the theoretical results appreciably. 

The author's analyses dealing with the effect of viscosity on the 
flow and the energy distribution also require some amplification. 
First of all, 7; and r; in Equations (72) through (77) are not de- 
fined. Are they the same as R, and R, in Fig. 1, referring to the 
inner and outer radii of the injection slots, or are they arbitrary 
radii? In the section labeled “Effect of Viscosity,’’ no symbol de- 
noting viscosity appears, but in Equation (65) € is defined as the 
turbulent exchange rate. This terminology is ambiguous, in 
particular since the equation for the shearing stress in turbulent 
flow is usually written 


dq _4 
dr r 
where u is the viscosity of the fluid and ¢€ is the eddy diffusivity. 

If the symbol ¢ in Equation (65) is intended to denote a com- 
bined coefficient of viscous and turbulent shear, it is difficult to 
reconcile this with the statement following Equation (77). ‘““The 
change from one flow to the other is due to viscosity, though it 
seems surprising that the relatively low viscosity of the air may 
have such an effect.’’ The validity of the analysis becomes even 
more questionable when one attempts to interpret the conclusion 
“that the energy of the rotational flow is equal to that of the irro- 
tational flow minus the energy taken up by the viscous forces.” 
Since the system used in the analysis, ie., the vortex tube, is 
adiabatic, no net transfer of heat into or out of the system can take 
place. Consequently, viscous or dissipative forces can produce 
changes in the energy and velocity distribution, but cannot alter 
the total energy content. If the word energy in the foregoing 
quotation is intended to denote the sum of the enthalpy and the 
kinetic energy stored in the fluid, one is at a loss to explain what is 
meant by “the energy taken up by the viscous forces.” The 
cross-sectional distribution of the energy stored in the fluid can 
only be evaluated when the static temperature as well as the ve- 
locity of the fluid are known. However, the static temperature 
distribution cannot be obtained from a macroscopic analysis of the 
flow, but requires a detailed analysis of the interaction between 
the radially outward transfer of energy as work and the radially 
inward transfer of heat. In the absence of such an analysis it is 
difficult to see how the author arrived at Fig. 10. 

Clarification of the questions raised in the foregoing discussion 
will be appreciated by those interested in applying the author’s 
experimental data and his theoretical analysis of the flow pattern 
to future studies of the problems associated with vortex-flow 
temperature separation, 


Author's Closure 


The author thanks Messrs. Gambill and Kreith for their writ- 
ten discussion. He also would like to express his appreciation to 
those persons who have communicated with him informally. 

Concerning Mr. Gambill’s discussion, the author acknowledges 
that »?/r? in Equation (3) should read v?/r, and the K2/r in the 
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equation following Equation (7) should read K?/R. As for the 
axial velocity distribution, the author assumed a constant axial 
velocity in Part 2 of his paper. This is somewhat the same 
assumption that Messrs. Diessler and Perlmutter made in their 
study of energy separation in vortex flows. There is at present 
very little detail pertaining to axial velocity distribution because 
of the difficulty of making accurate measurements inside a vortex 
tube. Part 1 of the paper reports a little on the distribution, and 
also gives a graph, but the author was not completely sure of the 
data, and he did not make use of it in the theoretical analysis in 
Part 2. The assumption of constant axial velocity is not en- 
tirely correct, but in the absence of detailed and reliable experi- 
mental data, and also for the sake of obtaining a mathematical 
solution, it was used in Part 2. There is at present a program 
being conducted at Michigan State for the Army Ordnance with 
the objective of getting more internal data. It is expected that 
in the near future, a more accurate analysis will be possible. 

Concerning Professor Kreith’s discussion, the author agrees 
with him that the assumptions made in the course of the mathe- 
matical developments are quite restrictive, and that the final 
results are likely to deviate from actual conditions. It must be 
remembered, however, that this is a first paper in the sense that 
it starts back from scratch and does not aim to take in too much 
territory. It is rather difficult to even get started on any analy- 
sis without some simplifying assumption. The model postulated 
in Part 2 is qualitatively verified by the flow pattern obtained in 
Part 1 (Figs. 20 and 23), at least upon entrance into the vortex 
tube. The quantitative verification by means of Equation (59) 
is a good suggestion. So far, it has not been attempted because 
it was felt that more data at different stations were desirable. 
The effect of viscosity was discussed in the last portion of Part 
2. The discussion, however, was not intended to be a full treat- 
ment of the subject. The author is aware of the shortcomings in 
the present treatment; he has thought of formulating a more 
accurate one, but so far, has not had much success. In Equa- 
tions (72) through (77), rm and r2 designate arbitrary radii in the 
vortex tube. They do not refer to the radii of the injection slots 
as Fig. | might have implied. Physically, rz would approach the 
radius of the vortex tube, and r; would approach zero. As for 
Equation (65), the author used Kassner’s notation for turbulent 
shear in circular flow, lumping the quantity 1 + pé€y into a single 
factor e. Admittedly, this is courting ambiguity, but it was done 
with the idea of reducing the number of terms to be carried about 
in future mathematical developments. It really doesn’t make 
much difference whether a combined coefficient of viscous and 
turbulent shear is used or not, since the eddy viscosity is so much 
larger than the molecular viscosity. The statement attributing 
to viscosity for the change from one flow to another is consistent 
with the theoretical developments in the paper when it is kept in 
mind that the viscosity referred to is primarily the eddy viscos- 
ity. As for the interpretation of the statement concerning the 
energy of the rotational flow and that of the irrotational flow, all 
that it means is a rearrangement in velocity and energy distri- 
bution. It does not imply an alteration in the total energy con- 
tent. It also aimed to point out the fact that irrotational flow 
(which exists at the start), by virtue of its more orderly arrange- 
mept of concentric streamlines, had greever “availability”? than 
rotational flow, which comes about through the mechanism of 
turbulent momentum exchange. 

It is hoped that these remarks wiil be of help in answering some 
of the questions pertaining to Parts 1 and 2 of the present paper. 
A better and more complete paper (based on experimental work 
now in progress) is being considered. In the meantime, the 
author thanks his correspondents for their expressions of interest. 
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a Re-Entrant Space Vehicle 


An analysis of the aerodynamic heating is presented through the trajectory and over the 
surface of a re-entrant hypersonic space vehicle. Bodies exhibiting zero and high lift 
over drag ratios are considered. The turbulent and laminar convective heat inputs are 


specified as functions of the trajectory and space vehicle parameters. The maximum 
heating rates and time integrated heat fluxes are given as functions of the local pressure 


distribution, body geometry, and wall temperature. 


Examples are presented to il- 


lustrate the application of this analysis. 


Introduction 


NE OF THE PRIME DESIGN considerations for a blunt 
nosed re-entrant space vehicle is the aerodynamic heating im- 
posed during re-entry. The heating characteristics such as the 
magnitude of the heat flux, total heating, and heating time will 
depend on the re-entry trajectory, the geometry of the vehicle, 
and the nature of the heating (turbulent or laminar). The heat- 
ing characteristics for a space vehicle will have a considerable in- 
fluence on the configuration and vehicle design. It is, therefore, 
necessary to establish an accurate evaluation of the aerodynamic 
heating so that proper considerations can be given to the heat 
protection system requirements for the space vehicle. 

As the vehicles of interest require hypersonic velocity during 
re-entry, the environment is one of extremely high temperature. 
During re-entry, therefore, compressibility and real gas effects 
must be considered in order to properly evaluate the magnitude 


of the aerodynamic heating imposed on the vehicle. It has been. 


demonstrated [1, 2]! that the reference enthalpy method can 
approximate real gas and compressibility effects. The reference 


! Numbers in brackets designate References at end of paper. 

Contributed by the Heat Transfer Division and presented at the 
Aviation Conference, Los Angeles, Calif., March 8-12, 1959, of Tue 
AMERICAN Society oF MECHANICAL ENGINEERS. 

Note: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those 
of the Society. Manuscript received at ASME Headquarters, Decem- 
ber 2, 1958. Paper No. 59—Av-5. 
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enthalpy method will, therefore, be used throughout the analysis. 

The aerodynamic heating analysis will be made for both turbu- 
lent and laminar-type boundary layers. The turbulent heating 
will be extended to all parts of the vehicle while laminar heat- 
ing will be directed primarily to the stagnation point. The prime 
purpose of this analysis will be to arrive at a systematic method 
for analyzing the re-entry heating for specific vehicle geometries 
in a form that can be applied readily to design. The analysis will 
be made primarily for three dimensional bodies, however, it can 
be readily extended to two dimensions. 

Since the trajectory has a major influence on the heating, two 
typical types will be used in the analysis. They will cover the 
majority of the cases of current interest and will also represent 
the two extremes in aerodynamic heating that may be encoun- 
tered. Zero and high lift vehicles will be considered. 

From this analysis a convenient and fairly accurate method of 
estimating the aerodynamic heating for a blunt nosed re-entrant 
space vehicle will be established which should give an understand- 
ing of the heating problems involved in space vehicle design. 
Previous papers [3, 4] have discussed the problem in a general way 
but this analysis will essentially present tools for the determina- 
tion of the aerodynamic heating of any given space vehicle ge- 
ometry. - 


Ballistic Vehicle 


1 Trajectory Considerations. The analysis of the motion of a 


drag over lift = D/L 

reference area, ft® 

velocity of sound, ft/see 

trajectory parameter, po/8(W/Cp 
A) sin O, 

wall temperature parameter, 
2gSh,,/U2 x, sat 


= drag coefficient 


lift coefficient 

specific heat, Btu/lb-deg F 

drag, lb 

acceleration of gravity, ft/sec? 

enthalpy, Btu/Ib 

mechanical equivalent of heat = 
778, ft lb/Btu 

lift, Ib 

mass, lb sec?/ft 

Mach number 

pressure, lb/ft? 

heat transfer rate, Btu/ft* sec 
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time-integrated heat-transfer rate 
(total heat), Btu/ft? 

distance of vehicle from center of 
earth, ft 

radius of earth, ft 

gas constant, ft/deg R 

nose radius, ft 

wetted length, ft 

time, sec 

temperature, deg R 

velocity, ft/see 

satellite velocity (reference), 
25,930 ft/sec 


= dimensionless velocity (U/Usat) 


weight, lb 

altitude, ft 

constant = 1/22000, ft 

ratio of specific heats (free 
stream) = 1.4 


ratio of specifie heats (edge of 
boundary layer) = 1.2 

local body angle, degree 

re-entry angle, degree 


= viscosity lb/ft sec 


partial range, radians 
density, lb/ft* 


total range, radians 


re-entry conditions 
burnout 

maximum 

reference conditions 
stagnation point 
turbulent 

wall conditions 


free stream conditions 
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hypersonic space vehicle with zero lift re-entering the earth’s 
atmosphere can be simplified by neglecting gravity without sub- 
stantially affecting the accuracy of the trajectory [3]. This is 
reasonable since the deceleration of the vehicle is large compared 
to the acceleration of gravity. The flight path angle will then be 
constant for this trajectory. The equation of motion for the 
vertical component may be written as: 


dty dU 


dl J 
dt dt 


dy 


sin O, 1) 
= ( 


2mg 


Based on the assumption of an isothermal atmosphere, the 
density may be written as: 


Pp = poe ~ By 2 
where 
po = 0.1095 lb/ft? B = 1/22000 


The velocity can then be obtained from equation (1) if Cy is 


assumed constant 


U = Ue 2 


where 


b = po/B(W/C,A) sin O, 


and the time from re-entry can be obtained as a function of altitude 
from equation (3) as 


1 E? 
+ 
BL’, sin O, 2-2! 3-3! 


where 


This trajectory information will be sufficient to make a detailed 
aerodynamic heating analysis. The variation of the velocity and 
altitude with re-entry time is shown on Fig. 1. The trajectory 
time and altitude are based on re-entry at 280,000 feet. The 
altitude, velocity, and time are shown as dimensionless quantities 
referred to the re-entry conditions of velocity and flight path 
angle. The curve shows the effect of the parameter b on the tra- 


jectory. The trajectory is shown for values of b of 1 to 50. 
12 T T T T T 
| VELOCITY | 
> | >" | | 
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~ 
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= 
i +- i iL 
2 4 6 10 20 22 


Ue SING, 
Fig. 1 Trajectory for zero lift vehicle showing dimensionless velocity and 
altitude versus time 


The lower b becomes, the longer is the time required before the 
velocity is materially changed from the initial re-entry velocity. 
At a given re-entry time, as b decreases, the altitude decreases and 
The trajectory parameter b 
is inversely proportional to W/C)A and sin Og as shown in equa- 
tion (3) 


the velocity increases significantly. 
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2 Turbulent Heat-Transfer Analysis. It should be noted at the 
outset that this analysis considers only the convective aerody- 
namic heating to the re-entry body. This will normally constitute 
the major portion of the heating, however, at hypersonic velocity 
the contribution due to radiation from the gas cap may become 
appreciable for some specific vehicles. 

The turbulent heat-transfer rate can be expressed conveniently, 
considering real gas effects, by a modified flat plate reference 
enthalpy method [1, 2] as: 

where 
j = (0.5h, + 0.28h + 0.22h,)/96.8 

This equation can be written more conveniently by expressing 
the quantities (p, U, h) behind the bow shock wave of a blunt 
nosed body in terms of free stream conditions (¢,, and LU’) by 
assuming a Prandtl number of unity. For a wall temperature of 
zero (h,./h, = 0) the turbulent heat flux can be written as: 

Ge = (6) 


where 


1 
f(MyP/P,) = 2.07 x 
0.4 


(P/P,) ] + 1.272(P/P,) 4 | 


For a given point on a body P/P, is constant at high Mach num- 
bers (M > 5), therefore, f(MyP/P,) is constant. 

The turbulent heat-transfer rate of equation (6) can be written 
in terms of the trajectory equations (2) and (3) as: 


3.48, —By 
— (0.88) be 


and the maximum heat-transfer rate is then 


Gem = S(MyP/P ~-*p04U (8) 
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Fig. 2 Maximum turbulent heat flux and local body angle versus local 
pressure ratio for a zero lift vehicle 


The maximum turbulent heat-transfer rate is shown in Fig. 2 as 
a function of the local pressure ratio. The heating rate function, 
F (4m), reaches a maximum at a pressure ratio of 0.6. For given 
re-entry conditions, body location, and local pressure ratio, the 
maximum turbulent heat flux can be readily evaluated. The 
pressure ratio is also shown as a function of local body angle for 
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Fig. 3 Maximum turbulent heat flux versus body position for a zero li 
vehicle 


assumptions of Newtonian and Newtonian plus Prandtl-Meyer 
theories [5]. It is noted that the Newtonian pressure distribution 
implies free stream pressure for a zero body angle. The zero in- 
tercept shown was based on a Mach number of 20. 

The heat-transfer rates can now be specified for any given body 
shape and re-entry condition. The maximum heat-transfer dis- 
tribution on the surface of a sphere cone has been evaluated for a 
one-foot nose radius with varied cone angles as shown on Fig. 3. 
For the purpose of this curve, the pressure on the conical portion 
of the sphere cone was assumed to be constant. The heating rates 
are given for both Newtonian and Newtonian plus Prandtl-Meyer 
pressure distribution assumptions. It is noted that the maximum 
turbulent heat-transfer rate occurs at a local body angle of ap- 
proximately 50 degrees. The dimensionless variation of turbulent 
heat-transfer rate to the maximum rate as a function of altitude 
can then be obtained from equations (7) and (8) as 


G:/Gtm = (2.175be -(0 ) (9) 


The time integrated heat-transfer rate can be conveniently ex- 
pressed in terms of the maximum heating rate as 


"sn. 
Q./din = dim dy y 


Ug sin 0, 


0 
Ug sin 0; Jo 


~Usin 0, 
Q1/dim = 


The integral in equation (10) can be expressed as an incomplete 
gamma function [6, 7], therefore: 


(1.752¢)4 


= ———— I'(038, 1.248 


(11) 


The time integrated total heat flux is then 


MyP 0.85 -0.2}, 


2) 
1.188 B sin 


It should be noted that the derived aerodynamic heating relation- 
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ships were based on assuming h,/h, = 0. When the wall tem- 
perature is greater than zero (h,,/h, > 0), the heat transferred to 
the body will be reduced as follows: 


— By 
2gJh,e* | 
Ith, =0 


1+ 12 


5-1 
72(P/P,) 
i=! 


4. 54gJh, 
L.272(P/P.) 


1+ 


A similar correction to the total aerodynamic heating can also 
be obtained by integrating equation (13) over the trajectory 
from re-entry to any desired altitude. 

3 Laminar Heot-Transfer Analysis. The stagnation point heat- 
transfer rate can be obtained from Lees’ analysis [8] modified to 
account for real gas and compressibility effects using the reference 
enthalpy method as follows: 


1.526 10-6 VU, 'P, Ah (14) 


The stagnation heating can now be written in terms of the re- 
entry conditions, trajectory parameters, body shape, and wall 
temperature as was previously described in the turbulent heating 
analysis. Equation (14) can be expressed in terms of free stream 
conditions [9] as: 


4, = f(My) ( ee.) (15) 
N 


where 


~ + 3\'/s 
f(My) = 6.035 x € (7; ') 
7 + 1 y¥+1 27 


The stagnation heat-transfer rate can then be written in terms 
of altitude as 


q, = f( My) (2 
N 


and therefore the maximum heat-transfer rate is 


Po 
3.2beR 


The dimensionless ratio of stagnation heating is then 


and the time-integrated stagnation-heat flux ratio is 


The integral in this equation may be expressed as an error func- 
tion, therefore: 


(16) 


Gum = S(MY) ( 


B 3.2 
be 
G./Gom = (3.2be)'*e G 2 


(3.2be)'/* 
Uy sin 6, 


Q, /Gum 


erf Vib 


4.98 
= 


(20) 
BU, sin Og 


The time-integrated heat flux is then 


sin 


T Po 


f 1.1bR, 


When the wall temperature is greater than zero, the stagnation 
point heating rate is reduced as follows: 
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The analysis given will establish the stagnation point heat 
transfer as a function of the trajectory parameters. In order to 
evaluate the laminar heat transfer at positions on the body, other 
than at the stagnation point, the following relationship of Lees 
[8] applies: 

4 1/2(P/P,(U/U.)R VR, 


=z 


0 | dO 


(23) 


The effect of shape on laminar heating is shown on Fig. 4 for 
sphere-cone bodies having various conical skirts. The pressure 
distributions of Newtonian and Newtonian plus Prandtl-Meyer 
are compared. The maximum laminar heating is always at the 
stagnation point on the body. 

4 Comparison of Turbulent and Laminar Heating. The variation of 
turbulent and laminar heat flux over the trajectory of typical 
zero lift vehicles is shown on Fig. 5. The dimensionless heating 
ratio (q¢/q,,) is given as a function of dimensionless re-entry time 
The maximum turbulent 
heating occurs at longer re-entry times, which implies lower alti- 


for various trajectory parameters (b). 


tudes and lower velocities, than that for maximum laminar heat- 
ing. This can be readily seen by examining equations (7) and 
(16) for maximum heating. As the trajectory parameter b de- 
creases, the heat flux curves are displaced to higher re-entry 
times. The magnitude of the laminar time integrated dimension- 
less heat flux curves (dashed lines on Fig. 5) appears larger than 
that for the corresponding turbulent case. This is clearly shown 
in Fig. 6 where Q/@,, is given for both the turbulent and laminar 
heating as a function b. The dimensionless total heating ratio 
(Q/q,,) is about 4.1 for turbulent and 5.0 for laminar heating when 
the trajectory parameter is greater than 5. The laminar dimen- 
sionless total heating ratio is, therefore, over 20 per cent greater 
than the turbulent ratio. 

The maximum heat-transfer rate correction for wall tempera- 
ture is shown on Fig. 7 as a function of the wall temperature 
parameter b. The effect of wall temperature on turbulent heating 
is greater than that for laminar. The correction for total turbu- 
lent and laminar heating is shown in Fig. 8 as a function of the 
trajectory parameter for various wall temperatures. For high 
values of b the correction is constant for a given wall tempera- 
ture. 


226 / aucust 1959 


1B Ue SIN 


Fig. 5 Dimensionless turbulent and laminar heat flux ratio versus time 
for zero lift vehicles 
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Fig.6 Total heating ratio versus trajectory parameter for zero lift vehicles 
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Fig. 7 Temperature correction for maximum turbulent and laminar heat- 
ing rates versus wall temperature parameter for zero lift vehicles 
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Fig. 8 Temperature correction for total turbulent and laminar heating 


versus trajectory parameter for zero lift vehicle. (Note: Solid line— 
turbulent; broken line—laminar.) 
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Glide Vehicle 


1 Trajectory Considerations. The motion of a hypersonic glide 
vehicle may be expressed as 
‘2 
L — mg cos 9 = = 
— D+ -mgsin 90 = mdU/dt (24) 


Based on the assumptions of small path angle, constant gravita- 
tional force, constant L/D, and U4. = gro, the equations of 
motion can be reduced [4] to: 

2aU? 


dU? _ |, 
ds To tind 


The free stream velocity then can be obtained as 
Ve =1-(1- 
where a = D/L 


The maximum range is then [4] 


1-V; 


and the burnout velocity is 
= 1 — (28) 


The free stream velocity can then be written in terms of the range 
parameter (a@) from equations (28) and (30) as 
= 1 — 0/9) (29) 


and the time of flight is then 


(30) 
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Fig. 9 Trajectory for glide vehicle showing dimensionless velocity and 
altitude versus time 


The variation of velocity and altitude with flight time is shown 
in Fig. 9. The dimensionless velocity shown applies to any con- 
figuration and is a function only of the range parameter a@. The 
altitude, however, is a function of W/C,A anda@. The altitude 
corresponding to a W/C,A of 500 was used as a reference. The 
burnout velocity (V,) will vary from 0.84 to 1.0 and the altitude 
frora 175,000 to 250,000 feet for range parameters of 0.5 to 2.0, 
respectively. The time is a function of the lift to drag ratio (L/D). 
The altitude required for flight increases with the range parameter 
but decreases with W/C,A. 

2 Turbulent Heat-Transfer Analysis. The anticipated environ- 
ment of a hypersonic glide vehicle is such that compressibility and 
real gas effects must be considered. The turbulent heat-transfer 
rate which has been derived previously in equation (6) can now be 
written in terms of the glide vehicle trajectory parameters from 
equations (24), (28), and (29). 
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2Wg 


j, = 


0.8 
[1 — (31) 


This then states the turbulent heat-transfer rate at any point of 
the trajectory for a given point on the body. The heat transfer at 
burnout is then: 


2Wg\°* ‘ 
diy = f(MyP/P,)S~ A ) — 206.04 
(32) 
The dimensionless heat flux can then be obtained as a function of 
range from equations (31) and (32) as 
4 — —0/6))0.94 
= = (33) 


(i - e~ 290 )0.94 


The maximum heat flux can be readily determined from equation 
(31) as 


dim = f(MyP/P,)S ( 
CLA 


3 U sar! | 0.778) (34) 


The time-integrated turbulent heat flux can now be obtained 


from 
0 


from equation (25) 


1 dU? 


2 ds 


= age~ 


and from equation (26) 


au —6/¢) 
10 =al sat e~ 2a9(1 9/%)) 


therefore substituting equations (36) and (37) in (35) 


U sat fila) 


= ag (1 — e~74%)0.94 


where 


0 


1.6 3.6 
= 0.125 fag ,5.6a¢ _ 
5.6 


The time-integrated heat flux is then 


2w 
Q, = fi MyP/P, ( U sat? filam (39) 
a 


Since these relationships were based on assuming the wall tem- 
perature zero, the following correction to the heat flux is necessary 


where h,/h, > 0 
1 + 1.272(P/P,) 7 


1 + 1.272(P/P,) 7 


wee + 2.27h,/h, 


(40) 


similar to the correction for the zero lift trajectory equation (13) 
except for a high lift trajectory the stagnation enthalpy cannot 
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conveniently be expressed in terms of altitude without introducing 
W/C,A. 

3 Laminar Heat-Transfer Analysis. The stagnation point heat flux 
can be expressed in terms of trajectory parameter from equations 
(28), (29), and (24) in (15) as 


2Wg 


[1 — (41) 


and at burnout it is 


C,ARy 


The dimensionless heat flux ratio is therefore 


Gey = S(My) ( 


The maximum stagnation heating is 


2Wg 


Usar™? — (44 
ant) 


dom = ( 
and the time-integrated heat flux ratio is then obtained from equa- 
tions (36), (37), and (43) as: 


Une 


where 


= f e**(1 — e ~ 240(1 0.6 dg 
0 


& — 1) — 0.2¢ 
— 0.024e~ 44% 1)... 


and then the total stagnation heating is therefore 


0.5 flag) 


The wall temperature correction is 


(47) 
dune (1 + hy/h, 

The stagnation point heating can now be readily determined for 
any point on the trajectory. In order to estimate the laminar 
heating at other points away from the stagnation point the method 
of equation (23) applies. 

4 Comparison of Turbulent and Laminar Heating. The aerody- 
namic heat flux comparison for the glide vehicle is shown on Fig. 
10. The dimensionless heat flux ratio (¢/¢,) is given as a func- 
tion of referred flight time for various values of the range parame- 
ter ad. The heat flux ratio is appreciably greater for turbulent 
than for laminar flow. The maximum turbulent ratio (4,,/q,) 
reaches 2.5 times the laminar ratio at a range (a@) of 2.0. The 
maximum turbulent heat flux occurs at higher flight times than 
laminar which implies lower altitudes and lower velocities. 

The burnout heat flux for both the turbulent and laminar case 
is shown on Fig. 11. The maximum turbulent heating occurs at 
a local pressure ratio of 0.6 and also at a range (a@) of approxi- 
mately 0.4. The maximum laminar heat flux occurs at a range 
(a) of about 0.6. At satellite re-entry (Vf = 0.99, a@ = 2.0), 
the heating is about '/, and '/; of the maximum for turbulent and 
laminar flow, respectively. The dimensionless total heating ratio 
is shown in Fig. 12 as a function of a@. The total! turbulent 
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Fig. 10 Dimensionless turbulent and laminar heat flux rates versus time 
for a glide vehicle 
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Fig. 11 Turbulent and laminar heat-transfer rates at burnout versus 
trajectory parameier for a glide vehicle : 
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Fig. 12 Total turbulent and laminar heating ratio versus trajectory 
parameter for a glide vehicle 
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Fig. 13 Temperature correction for turbulent and laminar heating rates 
versus wall temperature parameter for a glide vehicle 


LAMINAR 


O75 1.0 125 175 20 
TRAJECTORY PARAMETER ~ o@ 


225 


Fig. 14 Temperature correction for total heating versus trajectory 
parameter for a glide vehicle 


heating ratio is greater than the laminar and independent of the 
local pressure ratio. 

The correction for wall temperature is shown on Fig. 13 for the 
maximum turbulent and laminar heat flux and in Fig. 14 for the 
total heat. In a manner similar to the zero lift trajectory, the 
heat flux temperature correction is a function only of the wall 
temperature while the total heating correction is a function both 
of wall temperature and range. The temperature correction is 
greater for the laminar heating. The wall temperature has a large 
influence on the heat transfer to a space vehicle and must be con- 
sidered for proper evaluation of a heat protection system design. 


Conclusions 


The effect of the trajectory on the aerodynamic heating of a 
blunt nosed re-entrant vehicle is summarized in Table 1. There 
is a direct similarity between the heating for the zero lift and high 
lift bodies for all the parameters except the reference velocity. 
The aerodynamic heating varies as (W/Cp,,A) to the 0.8 and 
0.5 power and a physical dimension (S, Ry) to the —0.2 and —0.5 
power for turbulent and laminar flow, respectively. The heat- 
transfer rates vary roughly as the cube and the total heating 
varies as the square of the reference (re-entry) velocity for the 
zero lift bodies as shown in Table 1. For the glide vehicle the in- 
fluence of the re-entry (burnout) velocity on the heat flux and 
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Table 1 Comparison of turbulent and laminar heating 
Ww 


-——Laminar—-~ 
stagnation point 


Ex- 


Quantity ponent —~—Turbulent—— 


qm 


y 


00 00 


i? 


orbs 


OO 
once 
orm o 


A zero lift vehicle 
B high lift (glide) vehicle 


total heating cannot be as simply described. The maximum heat 
flux is practically independent of the burnout velocity [equations 
(34) and (44)] when the range parameter (a@) is greater than 0.39 
for turbulent and 0.58 for laminar flow. For lower values of a@ 
the maximum heat flux is a function of the burnout velocity since 
the heat flux at burnout is a maximum. For a given range of 
flight the total heating is a function of the re-entry velocity and 
the maximum total heating occurs when the burnout velocity 
ratio equals 0.75 and 0.89 for turbulent and laminar flow, respec- 
tively. It can be stated, however, that the heat-transfer rates 
are roughly proportional to the square and the total heating as 
the cube of the glide vehicle reference (satellite) velocity. 

In general, the two types of vehicles analyzed showed both ex- 
tremes in the aerodynamic heating environment. The zero lift 
vehicle is exposed to very short time, high heat fluxes, and 
moderate total heating during re-entry. The high lift vehicle on 
the other hand is exposed to an environment of long time, low 
heat flux, but very high total heat inputs. 

The analysis described in this paper evaluates the turbulent 
and laminar aerodynamic heating through the trajectory and 
over the surface of a re-entrant space vehicle for both zero and 
high lift bodies. It should prove useful for preliminary design 
purposes in evaluating the effect of various re-entry vehicle 
shapes on the aerodynamic heating. 
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A Photographic Study of Boiling 
in the Absence of Gravity 


A photographic study was made to determine the qualitative effect of zero gravity on the 
mechanism of boiling heat transfer. 
for boiling water and a high-speed motion-picture camera. 


The experimental equipment included a container 
To eliminate the influence 


of gravity, these were mounted on a platform which was allowed to fall freely approxt- 


mately 8 ft 


During the free fall, photographs were taken of boiling from various surface 
configurations such as electrically heated horizontal and vertical ribbons. 
was varied to produce conditions from moderate nucleate boiling to burnout. 


The heat flux 
The re- 


sults indicate that gravity plays a considerable role in the boiling process, especially in 


s DESIGN OF POWER PLANTS and heat-exchange 
equipment for satellite and space vehicles has stimulated an in- 
terest in heat-transfer phenomena under the conditions of a very 
small or, in many cases, negligible gravity field. Heat-transfer 
processes such as free convection and film condensation would be 
expected to be strongly dependent on the body forces due to the 
gravity field. The boiling process is especially interesting since 
the dynamics of bubble motion are quite complex and are proba- 
bly influenced by a number of factors. These include the buoy- 
ancy forces arising from the gravity field, and the forces re- 
sulting from surface tension and from the expansion of vapor into 
the bubble cavity. It is the purpose of this work to gain some 
insight into the importance of gravity in the boiling mechanism 
by making a photographie study of the boiling process in the ab- 
sence of gravity. 

The apparatus includes a high-speed motion-picture camera 
and a heated container for boiling water. These are both 
mounted on a platform which is allowed to fall freely a distance 
of about 8 ft, and motion pictures of the boiling phenomena are 
taken during this time. During the free fall, the inertia forces of 
the system exactly balance the forces due to the weight of the 
body, and hence the boiling takes place in a condition of apparent 
weightlessness. A variety of boiling conditions were observed. 
These include nucleate pool boiling from the bottom of a con- 
tainer, and boiling from both horizontal and vertical electrically 
heated ribbons with various heat fluxes including burnout, Only 
photographic information was obtained, and no attempt was 
made for this series of tests to determine numerical magnitudes 
of heat fluxes or heat-transfer coefficients. 


Experimental Equipment 


A schematic representation of the apparatus is shown in Fig. 1. 
The main component in the structure is a platform in the form of 
a cross, which supports a high-speed Fastex motion-picture cam- 
era and a container for boiling water. At the ends of three 
branches of the platform there are large guide holes through which 
'/-in. vertical cables pass. The cables serve as a guide for the 
falling platform in case of an accident, and during a normal experi- 
mental run little or no contact is made between the platform 
and the cables. A balancing yoke is attached at the top of the 
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connection with the motion of vapor within the liquid. 


platform and can be adjusted so that the solenoid grapple, which 
holds the platform from above, will be located above the center 
of gravity of the freely falling equipment. 

The features of the apparatus are best described by consid- 


ering a typical experimental run. A 400-ml beaker is carefully 
cleaned and partially filled with distilled water. It is then 
placed on a hot plate which consists of a copper block mounted 
on a 500-watt electric heater controlled by a Variac. If boiling 
is to be studied from a horizontal or vertical ribbon the elec- 
trical terminals for the ribbon are fastened to the side of the 
beaker and the ribbon mounted at the desired location in the 
water. The ribbon power was supplied by a separate 18-amp, 
115-volt Variae. The ribbons used, were Nichrome, 0.006 in. 
thick, 0.125 in. wide, and approximately 7/s in. long. Power was 
then supplied to both the ribbon and the hot plate under the 
beaker, and the water was allowed to boil to remove absorbed 
gases. During the runs where boiling was observed from the 
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Schematic diagram of apparatus 


Transactions of the ASME 


pals 
Introduction 
7 
/ 
r 
3 
‘ 
| 


TIME, SEC 


Fig. 2. Pool boiling from the bottom of a glass beaker 


ribbon, the heater under the beaker was only used to help main- 
tain the liquid near the saturation temperature and the power 
was such that very little bubble formation occurred at the bot- 
tom of the vessel. The platform was then hoisted eight feet, and 
final adjustments were made on the heater power to obtain the 
desired boiling condition. 

The drive motor on the camera and the solenoid release grap- 
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ple which supports the platform are connected to a 6-sec indus- 
trial timer. The timer will start the camera and then activate 
the solenoid release in sequence. Generally, the camera was 
started 0.3 sec before the platform began to drop, and this gave a 
record of the initial boiling state before the gravity field was re- 
moved. The camera speed was generally 3500 to 4000 frames 
per sec, which gave a total running time for 100 ft of film of about 
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Fig. 3 Boiling from horizontal electrically heated ribbon with moderate heat flux; heat applied at beginning of free fall 


1.25 see; 0.3 sec initially, 0.7 sec during the free fall, and 0.25 sec 
during the deceleration at the end of the run. During some runs 
the boiling from the ribbon was started at the beginning of the 
free fall to eliminate any influence of a previous boiling process 
under normal gravity conditions. This was accomplished very 
simply by connecting the ribbon heater into the solenoid circuit. 

The platform was decelerated by vertical !/:-in. pipes, fastened 
to the bottom of the platform, plunging into a bed of sand 3 ft 
deep. The sand bed was found to be quite effective in bringing 
the apparatus to a gradual stop without damage to the camera. 
The camera was, of course, well supported on the platform. 

Adjacent to the boiling container a vertical measuring scale was 
positioned, so that a continuous record of the vertical location of 
the platform was recorded on the film. Also the camera re- 
corded timing marks on the film at a rate of 120 per sec so that a 
record of distance versus time was obtained. This indicated 
that the friction in the system was negligible and that a true free 
fall was obtained. 


Discussion of Results 


The results are obtained as motion pictures and it is naturally 
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more instructive to view these directly. A loan copy of the film 
is available on request from the authors. For the purposes of 
this report sequences of pictures were selected from the film and 
will serve to illustrate fairly well the description of the boiling 
process. 

Fig. 2 shows a series of pictures for pool boiling from the bottom 
of a glass beaker. The first two pictures in the sequence illus- 
trate the normal boiling process, and the third picture is at the 
beginning of the free fall. From observing the motion-picture 
films it was noted that there is on'y a slight inertia associated 
with the motion of bubbles through the liquid. Following the 
beginning of the free fall, the bubbles still continued to move a 
short distance upward and then remained suspended in the liquid. 
From the sequence of pictures it is seen that large bubbles are 
formed and increasing quantities of vapor are contained within 
the liquid during the zero gravity process. The liquid surface is 
less stable, and as evidenced by the last few pictures in the se- 
quence, a liquid film has been pushed up the side of the beaker. 
It was thought that this film might be in part due to surface ten- 
sion forces at the liquid glass interface which were no longer bal- 
anced by the gravity force. To examine this, a few runs were 
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Fig. 4 Boiling from horizontal electrically heated ribbon with high heat flux and burnout 


made using an unheated container, and during the 0.7-see dura- 
tion of the free fall, the liquid did not move appreciably up the 
side of the container. However, during boiling the agitation due 
to the formation and growth of the bubbles pushes liquid up the 
side and it remains there until the platform is decelerated. 

Fig. 3 shows boiling from a horizontal electrically heated ribbon 
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with a moderate heat flux. To eliminate the influence of bubbles 
formed prior to the free fall, the heating was started at the same 
moment that the platform was released, so all bubbles were 
formed under zero-gravity conditions. The hot plate under the 
beaker was used to maintain the fluid near the saturation tem- 
perature and a few bubbles were formed at the bottom of the con- 
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Fig. 5 Boiling from vertical electrically heated ribbon with high heat flux 


tainer, but these did not interfere with the vapor formation on the 
ribbon. The edge of the ribbon is viewed in the pictures, and 
the flat heating surface extends '/s in. into the plane of the paper. 
Under zero-gravity condit’ ms the bubbles remain attached to the 
ribbon and grow fairly large so that large areas of the ribbon come 
into contact with vapor. The size which the bubbles attain 
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during the duration of the fall is of course dependent on the mag- 
nitude of the applied heat flux. The bubble formation is 
equally extensive from the bottom surface of the ribbon as 
from the top. 

In Fig. 4 we again have boiling from a horizontal ribbon, but in 
this case the heat flux has been raised so that the ribbon burns 
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Fig. 6 Boiling from vertical electrically heated ribbon with high heat flux and burnout; heat applied at beginning of free fall 


out during the free fall. Also the ribbon is heated prior to the 
drop, and hence, the first few pictures illustrate the normal 
boiling process under high heat-flux conditions. During zero 
gravity a single very large bubble is formed roughly 1.5 in. in 
diameter (beaker diameter is 3 in.). The bubble in this case is 
located above the ribbon and this may be due to the initial 
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boiling process which caused the water above the ribbon to be at 
a higher temperature than the water below it. The last picture 
in the sequence shows the ribbon after deceleration. The large 
bubble which obscured the view of the ribbon has moved out of 
the liquid and the burned out ribbon can be observed. 


Figs. 5 and 6 show boiling from a vertical ribbon. The posi- 
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tion of the ribbon between the terminal posts is seen in the first 
picture of Fig. 6. In Fig. 5 the heat flux is at a high level, about 
10 to 20 per cent less than that required for burnout. The re- 
moval of gravity causes a single large bubble to be formed which 
is about 2 in. in size. e 

In Fig. 6 the heat flux has been raised so that the ribbon burns 
out during the free fall. The heating process is initiated at the 
beginning of the drop. At the right side of the beaker is a col- 
umn of bubbles originating from a nucleation center at the bot- 
tom of the container, and these bubbles remain suspended in the 
fluid during the run. In the picture at 0.633 sec the lower half 
of the ribbon can be seen glowing, at the beginning of burnout, as 
a slightly curved white line. 


Conclusion 


The photographic study indicates that the gravity field plays 
a very important role in the removal of vapor from the surface in 
a normal pool type boiling process. In the situations observed 
here, the fluid was always near the saturation point and the heat 
fluxes were sufficiently high to form large quantities of vapor. 
During the free fall the vapor remained adjacent to the heating 
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surface, and there was no evidence of bubbles being pushed away 
from the surface to any appreciable extent during their formation. 
The accumulation of vapor near the surface should make the 
electrical heater more sensitive to burnout, but the information 
obtained in this work is insufficient to make any quantitative 
statement on this hypothesis. In a situation where the liquid is 
highly subcooled, so that under ordinary conditions the bubbles 
form and collapse very close to the heated surface, the removal of 
gravity would be expected to have much less effect. 

Under some circumstances the process with zero gravity may 
actually be a transient one. If the liquid is at the saturation 
temperature or slightly superheated, the bubbles at the heating 
surface may continue to grow indefinitely, and hence larger bub- 
bles would be expected in an experiment with longer free fall 
times. If the liquid is slightly subcooled there is a chance that 
a steady-state equilibrium can be reached. Under these cireum- 
stances the bubbles would increase in size until their surface area 
is large enough to transfer to the surrounding subcooled liquid the 
heat being added by the heated surface. Under these conditions 
the bubbles would remain at a constant size. It appears that a 
situation like this was reached in Fig. 3. 
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= REJECTION OF SURPLUS HEAT during space flight 
requires possibly for the first time the design of radiators that are 
strictly radiators. 

Radiator configurations for space vehicles may be classified as 
single-surface or double-surface. The former type would be built 
into the skin of the vehicle (Fig. 1); the latter would be inde- 
pendent assemblies carried in the plane of the vehicle axis (Fig. 
2). In both cases the vehicle should be oriented (at least during 
periods of heat rejection) to keep the radiator surface in the same 
plane as the solar vector, unless shielded from solar radiation by 
the vehicle. Double-surface radiators might be nose-mounted, as 
indicated in Fig. 2(a), position 1, for minimum energy intercep- 
tion by the vehicle. 

Single-surface radiators, although economical in terms of strue- 
tural mass, suffer in effectiveness because of less ‘‘active’’ surface 
than two-sided designs; also they create insulation problems at 
the inner surface. On the other hand, double-surface radiators 
are likely to contain concavities or interior corners occasioned 
by minimum mass designs. Fin portions, material connecting 
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1 Study performed at Rocketdyne Division, North American Avia- 
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Nomenclature 


Some Radiator Design Criteria for 
Space Vehicles 


Nonconvective radiators of single or double-active surface design are analyzed for surface 
temperature distribution; also for view factors in the case of fin-and-tube geometry. 
Methods and examples of maximizing heat rejection per unit weight are given. 


tubes, Fig. 2(6), of these radiators have lower surface tempera- 
tures than metal in contact with the heat-transfer fluid; indeed 
the fin-temperature distributions are significantly influenced by 
view factors, as illustrated in Appendix 3 (Fig. 12). An analytic 
evaluation of geometric view factors of tubes and fins, for design 
of Fig. 2(), is made in Appendixes 1 and 2, respectively. 


Radiator Analysis 


From an analytical point of view three general types of radia- 
tors have been considered in this study. These are: 


1 Radiators whose surfaces have a constant temperature in 
the direction normal to the direction of the flow of fluid (07'/dz 
= 0) and whose geometric view factor F is 1.0, Fig. 1(b). 

2 Radiators which have intermittent fluid passages and, there- 
fore, cyclic variation in the surface temperature in the direction 
normal to that of fluid flow, and whose geometric view factor is 
1.0, Fig. 1(e). 

3 Radiators which have intermittent fluid passages and whose 
view factor is a function of position on the surface, Fig. 2(b). 


1 Flat Parallel-Plane Radiators (Radiation From One or Both Sides)- 
F = 1.0 and OT/0x = 0. Consider thermal radiation to an en- 
vironment of zero degrees Rankine from a radiator consisting of 
two parallel planes enclosing a flowing fluid, Fig. 3. Incident 
solar radiation has been assumed zero by orientation. Re-radia- 
tion and the interception of radiation by the supporting vehicle 
have been ignored. 


surface area 


fin width 


‘2 
0 V1 — sin? @ sin? 6 


radiator weight 
co-ordinate in radiator plane 


constant of integration 

specific heat at constant pressure 
of radiator fluid 

outer diameter of tube 

dimension of radiator normal to 
flow 

fin view factor 

view factor in plane normal to 
tubes 

tube view factor 

view factor at a position @ on 
the tube 

convective heat-transfer film co- 
efficient 

heat of vaporization of fluid 

thermal conductivity of radiator 
material 


oe/kt 
ce/kt 
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dimension parallel to direction 
of flow 

normal to tube or fin 

2(1 + B/D) 

heat rejection rate 

fin thickness 

tangents to tube surface 

T,, = f(z, y), surface tempera- 
ture 

tube surface temperature at a 
position y 

mid-fin surface temperature at a 
position y 

fluid bulk temperature 

fourth root-mean-fourth power 
tube surface temperature, 
equation (9d) 

Ty Tr 

dT /dr 

a?T 

fluid flow rate per tube, lb/sec 


normal to direction of flow 

co-ordinate parallel to direction 
of fluid flow 

co-ordinate normal to plane of 
radiator 

x — D/2 
B/2 

angular position on tube surface 

total emissivity of radiator sur- 
face 

fin effectiveness for thermal ra- 
diation 

angle between plane normal to 
fluid flow and other planes 
containing normal line N 

Boltzmann constant 

T/T xy 

angle of obstructed view from fin 
or tube surface 

— = view angle from 
normal to surface 
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C) POSSIBLE DESIGN 
Fig. 1 Single-surface space radiators 


FLUID FLOW 


From a steady-state heat balance on the elemental volume of 
width dy, Fig. 3, the heat radiated from the surface equals the 
heat lost by the fluid. 


—oeT'Edy = AT, (1) 


T, = 
h, 


If the radiator is constructed of a thin material of high thermal 
conductivity and if the convective heat-transfer film coefficient on 
the fluid side is large, i.e., if 


where 


t 
+ oer + T 


then 


aT, dT (2 
Za) 
dy dy 


Substituting (2a) into (1) and integrating between y = 0 and 


y=L: 
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Fig. 2 Double-surface space radiators 


(3) 


If the fluid ix condensing and the surface temperature can be 
assumed constant then 


wH , 
oeT' 
where H,, is the heat of evaporation. 


The total heat rejected is, Q = v@,(Tw — Ty) (4 


(or wH,,). 

These equations (2), (3), and (4) are sufficient to estimate and 
minimize sizes and weights of radiators of the first classification. 

2 Intermittent Passage, Flat or Convex-Surface Radiators (Radiation 
From One Side). F = 1.0and OT/0x = f(x). Consider thermal radia- 
tion from the flat (or convex) side of the radiator geometry shown 
in Fig. 4, assuming no solar or other incident radiation. Assume 
that heat transferred to or from the opposite side is negligible. 
Assume also that 07'/Oz = O from z = 0 tox = D/2 and therefore 
T = f(y) = forO < < D/2 and 0 < v < L, and that T = 


B+D 
f(z, y) = T,,for andO<y< L. 
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From a steady-state heat balance on the surface area element 


) dy the heat radiated equals the heat lost by the fluid 


flowing in the passage: 


D 


B 


where 


B+D 


f 2 
D/2 


B/2 


(6) 
and it is assumed that an appropriately determined value of 7 is 
a constant, 7, over the interval of integration in the y-direction. 


Making the substitution of equation (1a) and integrating from 
y= 0, T = Ty to y= L, T = Ti1, 


if the fluid is condensing and it can be assumed that 


aT, dT, _ 


dy dy 
then 


wH 
= fo (7a) 
+ 7B) 
It is now necessary to evaluate 9. Consider a steady-state heat 
balance on the volumetric element (dzrdy in Fig. 4. If we assume 
that conduction in only the z-direction is important, then 


—k dy = —k 
dx dr 


dT 
T + i - dz) dy + o€T*drdy 
ar 


dx? 


(8) 


Integrating and applying the boundary conditions (. r= 


B+D IT 
T = T>,, and 0), 
2 dx 


dT laa % 
= 404 (TS — 7:,'). (8a) 
dr kt 
By substituting with equation (8a) for dr in equation (6) and 
integrating equation (6) from T = 7, to T = To,: 


6 ye — 


oe 


and kK, = 72,3 


nB = 


If we let 


then equation (9) becomes 


oe 


ty 


nB = 2 
It will be shown later that 
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In other words, the fin effectiveness for the radiator is evaluated 
at the fourth root of the mean fourth power of 71,. 

An expression must now be derived relating the fin half-width, 
B/2, to the fin temperature drop, T7;, — T2,. Making the sub- 
stitutions of equation (9a) and rearranging, equation (8) be- 
comes 


ltor = 


B dt 
K,— = 
f Vr 


Equation (11) was evaluated graphically and the resulting plot of 
K,B/2 versus 7; is shown in Figs. 5 and 6. Therefore, given k, ¢, 
and ¢€, knowing 7), and 7.2,, B/2 can be found directly. When 7, 
and B/2 are known, 72, can be found by a trial and error method. 
Using equations (7), (9c), (11), (2), and (4) the size and weights 
of radiators of the second classification can be estimated and 
the weights can be minimized. 

The recommended procedure for using the previously derived 
equations for the design of a minimum weight radiator can best 
be described by working a typical example. 

Consider a radiator of the type shown in Fig. 4 which is cooling 
0.1 lb of hydrogen per second, per radiator tube from 2000 to 
1500 deg R (¢,, = 3.4 Btu/Ib deg R). Assume that structural and 
fabrication requirements necessitate a fin thickness, ¢, of 0.015 in. 
and that D = 0.5in., € = 0.95, k = 30 Btu/hr ft deg R, density = 
0.285 lb/in’ and the convective heat-transfer film coefficient on 
the fluid side, h, is 2000 Btu/hr ft? deg R. Determine the di- 
mensions, B and L, which will reject this amount of heat with the 
minimum weight of radiator structure (neglect the fact that part 
of the weight might ordinarily be considered vehicle structure 
weight). Note that in this example the material modulus, kt, is 
specified; however, in general kt should also be allowed to vary 
within the structurally prescribed limits when determining the 
minimum radiator weight. 

It has been found that by determining a single value of 7 over 
the length, L, evaluated at 7),,, gives sufficient accuracy and con- 
servative results as compared with an incremental method. From 
equation (2), Ti = 1986 and 7), = 1496 deg R, and therefore 
by equation (9d), Ti, = 1789 deg R. 

If the fin width B varies such that 7,, — 72, takes on values of 
50, 100, 200, and 300 deg R, the corresponding values of 7B can 
be determined using equation (9). Then from equation (7) and 
Fig. 5 or 6 we can determine L and B, respectively. Thus the 
radiator is geometrically specified and for each value of 7), — 
Ts the radiator weight per tube can be calculated (Table 1). 
By plotting radiator weight and B/2 versus 7:,, the minimum 
weight and corresponding fin width which will cool the hydrogen 
can be determined (Fig. 7). 

From Fig. 7 the minimum weight of radiator is 5.035 lb/tube 
where B/2 is 0.31 in. (B/D = 1.24) and = 1656 deg R. Note 
that this does not correspond to minimum surface area. 

Greater accuracy can be obtained if the radiator is analyzed in 
longitudinal increments of equal heat rejection. If ¢, is assumed 


B 
Integrating from z = — 
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Fig. 5 K,8/2 of K.B/2 versus + 


constant over the temperature range this then becomes incre- 
ments of equal fluid bulk temperature drop. This was done for 
the case of Table | where B/2 = 0.410 in. Each of these inere- 
ments was handled in a manner similar to the example given ex- 
cept that 7,,, was found from Fig. 5 or 6 by a trial and error proc- 
ess for each increment. Knowing 7;,, and B/2, a value of T.,, is 
assumed, thus enabling the computation of K,B/2. From Fig. 5 
or 6, knowing 7, A, B/2 can be read off the ordinate and checked 
with the previously computed value. This process is repeated 
until satisfactory agreement has been reached. The increment of 
equal fluid bulk temperature chosen was 100 R and the results 
are shown in Table 2 

In comparing the results of the iterative analysis, L = 46.26 ft, 
with that determined in the example of Table 1 (B/2 = 0.410in.), 
L = 47.0 ft, it can be seen that the latter is slightly larger and 
therefore represents a conservative estimate. Based upon this 
comparison it is recommended that in general the former method 
be used. 

Double-Surface Radiators: Fin-and-Tube Construction. fin- 
and-tube configuration lends itself well to a double-surface, high- 


(500 
Fig. 7 Example of radiator optimization 


three-dimensional effects, have been computed for this design for 
a wide range of the ratio B/D (Appendixes 1 and 2). A design, 
Fig. 2(c), approximating the fin-and-tube arrangement, which 
might be easily fabricated and satisfactory for low or medium 
pressure fluids, would have higher view factors. Hence analyses 
on the configuration of Fig. 2(b) would be conservative for this 
design 

In this analysis it is assumed that at a given position in the y- 


pressure, minimum-weight 


T, 
deg R 
1789 
1780 
1780 
1789 


deg R 
2000 
1700 
1600 
1500 
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design. 


T 2m, 
deg R 
1739 
1689 
1580 
1480 


nB, 
it 
0 0204 
0.04)0 
0.0541 
0 0629 


View factors, considering 


L, 
ft 
63.3 
54.4 
417.0 
43.1 


Weight. 


Ib/tube 


5 31 
5.06 
5.00 
5 34 


Heat rejected per tube = 17 Btu/sec = 17.92 kw 


Table 2 Results of iterative analysis (8/2 = 0.410 in.) 


Tw, 
deg R 
1791 
1693 
1594 
1496 


Tim: 


deg R 
1940 


Tm, 
deg R 


1690 
1630 
1560 
1495 
1416 


0B, 
ft 


0 0506 
0.0518 
0 0548 
0 0556 
0 O588 


AL, 
ft 


5.67 
7 
8 41 
10.80 
13.82 


direction the temperature on the surface of the tube is constant 
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ft? /tube 
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Fig. 8 Fin-and-tube construction 


(in the circumferential direction), i.e., on the tube surface T = 
f(y) = T,, and that the surface temperature on the fin is a function 
of xand y, T = T,,. 

From a steady-state heat balance on the volumetric element of 
width dy (Fig. 8), neglecting solar and other incident radiation, 


B+D 


2 


c,dT, = —oe F,T,, dz \dy, 


where F, is constant and F,, is a function of z (vp szs 


B+ D 
9 


) for a given B/D ratio. Integrating as in Section 2, 


— 2t)F, + 7B] = 


where 


7 is again evaluated at 7),,. 

As a result of a heat balance on the volumetric element tdrdy 
(Fig. 8) with thermal radiation from both sides, it can be shown 
that 

= 
dr? kt 


with the following boundary conditions: 


Tin (known 
(must specify 


and (symmetry at mid-fin) 


Because of the strong variation of F, with z, the solution of 
equation (13) requires a rather tedious trial and error step inte- 
gration (Appendix 3) and the result is a plot of F,7* versus 


z — D/2 
BP (for a typical example see Fig. 11). 


type 7 can be determined graphically. 
referring to Fig. 11 


From a plot of this 


From equation (13) and 
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Fig. 9 Fin temperature drop, T, — T. versus B/D (fin-and-tube construc- 
tion, condensing fivid) 


1000, 
T.=1800°R 


k= 30 


im 


“Yr 


Fig. 10 Heat radiated from fin per unit longitudinal length (fin-end-tube 
construction, condensing fluid) 


Fig. 11 Fin effectiveness versus position on fin as influenced by fin 
view factor (T; = 1800 deg R, T, = 1600 deg R) 
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(13a) 


Thus with equations (2), (4), (12a), (13a), and the methods of 
Appendix 3, radiator sizes and weights of the third classification 
can be estimated and the weights minimized. 

4 Double-Surface Radiators: Fin-and-Tube Construction. (Approxi- 
mate Solution Assuming F, = 1.0.) Because the methods of Ap- 
pendix 3 are so unwieldy and time consuming it was deemed neces- 
sary to present a simpler method for estimating sizes and weights 
of radiators of the third classification. By assuming that the fin 
view factor F, equals one, equation (13) reduces to the same 
easily solvable form as equation (8). In comparing the results of 
an example worked below using the methods of this section with 
the exact solution worked in Appendix 3 it is seen that assuming 
F, = 1.0 results in a 10 to 20 per cent overestimation of heat re- 
jections by the fin portion of the radiator, Fig. 10. It is proposed 
that, as a result of future analytical work, a generalized correction 
factor could be generated (possibly a function of B/D) which 
could readjust the fin heat rejection rate, as calculated by the 
approximate method, to show closer agreement with the exact 
solution. 

Assuming Ff, = 1.0 in equations (12a) and (13) and integrating 
equation (13) we find 


. 
7B = — 39 kt (Tra* — (14) 
equation (11) becomes 
ri 1 
where 
T= — and = yes (15a) 
Tx kt 


and, therefore, Figs. 5 and 6 are applicable. 

The use of equations (12), (14), and (15) for estimating the 
sizes and weights of radiators of the third classification can again 
best be demonstrated by means of an example. Consider a fluid 
(of negligible density) condensing in a fin and tube type of radia- 
tor under conditions such that the tube outer surface temperature 
is constant, 7,. Assume 7, = 1800 R, k = 30 Btu/hr ft R, both 
tube and fin thickness = 0.03 in., D = 0.5 in., and € = 0.95 
Equation (12) can be rearranged to give the heat rejected per foot 
of length, 


H 
Q/L = = = — 2)F, + 


4 


(16) 


Considering a section which includes half a tube (Fig. 8), the heat 
rejected from the tube section is 


aD 
(Q/L), = ( oeF ,T (16a) 
and the heat rejected from the half fin is 
nB 
(Q/L), = = ktKeT, — 1. (16b) 
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Fig. 12 Temperature distribution on fin (T, = 1800 deg R, T, = 1600 deg 
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Fig. 13 Heat rejection per Ib of radiator versus B/D (T,; = 1800 deg R, 
+ = 0.03 in.,D = 0.50in., k = 30 Btu/hr ft deg R) 


If we consider a range of fin temperature drops, 7, — 72, the 
corresponding fin widths, B/2, can be found using Figs. 5 or 6, 
the heat rejected per unit length can be computed from equation 
(16) and the weights can be calculated, knowing the geometry. 

The relationship obtained between (T; — T,) and B/2 is 
shown in Fig. 9 and is compared with the exact result from Ap- 
pendix 3. The heat radiated per unit length from a half fin is 
shown in Fig. 10, along with the corresponding result from the 
exact solution. The results are summarized in Table 3 and Fig. 13 
where the maximum heat rejection per unit weight is about 15,220 
Btu/hr, lb for B/D = 1.45 for the approximate solution and 
14,100 Btu/hr, lb for B/D = 1.55 for the exact solution 


Table 3 Fin-and-tube design summary (7, = 1800 deg R) 


Ts K:, B/2, 
deg R K.B/2, ft-! in. 
1750 1.029 0.135 9.64 0.168 
1700 1.059 0.196 9.25 0.254 
1600 1.125 0285 844° 0.405 
1500 1.200 0354 7.67 0.554 
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Q/L, Q/W, 

(= ) W/L, ) 
B/D hr ft lb/ft hr lb 
0.672 1320 0.0926 14,250 
1.016 1521 0.1016 15,080 
1.620 1779 0.1170 15,210 
2.216 2003 0.1323 15,140 
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Fig. 14 Determination of view factors on tube 


APPENDIX 1 
View Factors on the Tube 


In Fig. 14, the plane view of space from a point A on the tube 
surface is bounded by the tangents AT; and AR. The follow- 
ing expression gives the complementary angles W and @ (between 
AR and the normal, AN; and AR and the tangent AT», re- 
spectively) 


sin ¥ = cos ¢ = P sin @ V p? — 2p cos a — (p cos — 1) 
p? — 2pcosa + 1 
(17) 


where p = 2(1 + B/D) and a@ defines point A. To find the three- 
dimensional view factor we consider an angle 6 in the A-N plane 
(normal to the plane of the sketch) and measured from the 
plane of the sketch (Fig. 14). 


The space view factor for any 
position A is then 


cos @. 
K(®) is an elliptic integral 


dé 


K(é) « 
Vi — sin? @ sin? 6 


The quadrant, 0 < @ < 2/2, was divided into six equal sectors 
(a = 7.5 deg, 22.5 deg, ete.) and the average value of the view 
factor for the tube was determined and is shown in Fig. 15. It is 
interesting to note the comparatively narrow limits of tube view 
factor: at B/D = 0, F, = 0.648 and at B/D = o, F, = 0.873. 


(19) 


APPENDIX 2 
View Factors on the Fin 


The two-dimensional view of space from point z on the fin is 
bounded by tangents 27; and rT». The radiation received at a 
point S is proportional to cos Y = sin @. The view factor in the 
plane of the sketch (Fig. 16) is the ratio of the radiation emitted 
between the obstruction angles @; and @. to that emitted by an 
unobstructed plane 


r/2 
f sin ddd + f sin ddd 
2 f, sin ddd 
{cos + cos (20) 


To obtain the complete three-dimensional view factor we con- 
sider the view limits ¢ and @» established by intersections of 
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| | 
Fig. 15 Tube geometric view factor versus B/D (view of space) 
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Fig. 16 Determination of view factors on fin 


several planes (including the normal zN and making angles 6 with 
the plane of the sketch, Fig. 16) with the tangent planes x7’, and 
xT;. Integration of these view factors F's (expressions similar to 
equation 20) gives the complete fin view factor. The cut-off 
angles @; and @: in a plane defined by @ are determined by 


tan = tan d cos 0, tan de = tan gd cos@; (21) 


whence 
cos de = 
Vi + tan*® cos? 


1 
cos be = = . 
V1 + tan? cos? 


Replacing cos @ and cos ¢, in equation (20) by the expressions of 
22) and integrating on @ through a quadrant, 


= 
dé 4 2 dé 
0 V1 + tan? @; cos? 6 0 v1 + tan? od: cos? 8 


w/2 
2 
0 


(23) 
The denominator of equation (23) “normalizes” the view factor 
(amounting to evaluation of the numerator for “full view” @ = 
¢: = 0). Each term in the numerator of equation (23) may again 


be expressed in terms of a “complete” elliptic integral (ef, equation 
19): 


cos 2 dé 
0 V1 sin? sin? 6 


These substitutions yield 


K(,,2) 


1 
F, = cos + cos 
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Fig. 17 Fin view factor F, (at various positions along the fin X) versus 
B/D 


which has the expected extremes F (g, = @: = 0) = 1.0, P(g = 
= 4/2) = 0. 

The cut-off angles ¢, and q@» are functions of B/D and X = 
D/2 

B/2 


l 


sin p = i + X B/D’ 


1+(2—X)B/D 


sin @ = 26) 


The fin view factors are plotted as functions of B/D for various 
X in Fig. 17. It is interesting to note (Fig. 17) that as B/D— oa, 
only the view factor for X = 0 (next to the tube) approaches 0.5, 
while these at all other positions approach 1.0. 


APPENDIX 3 


Determination of Exact Temperature Distribution on Fin 
(for Fin-and-Tube Design) 
As brought out in Section 3, the temperatures on the fin are 
influenced by view factors (Appendix 2) which are functions of z. 
Repeating equation (13) and its boundary conditions, 


eT ) 
T°’ = =? F,T*; 

kt | 

at z= D/2, T = Tim (known), | 
> (13) 

at T = (specify), 

B+D dT | 

and = (by symmetry). | 

2 dr | 


One finds the somewhat novel situation here of knowing two 


+ D 


, at which 


boundary conditions at an unknown location, 


the step solution must nevertheless begin Since the purpose of 


these calculations is to provide a fin effectiveness or its equivalent, 
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Table 4 Exact solution for fin temperature distribution 


+B 0.00001 0.090 0.180] 0.270] 0.360} 0.450 


T (*R) 1600 1607 .€£ [1630.89 670.50 |1727. 1804.65 


(or)* 


Fig.17) 

2.72] 2.74 2.92] 9.16 3.40 2.28 
“AT (%) | 2077 | 2184 | 2368 | asse 


AT=-T dx CR) 


each end -QT = "dx, end - OT T’ elsewhere. 


they should be made at the fourth-root-mean fourth power of the 
tube surface temperature extremes, T),. 


Starting at —>—, mid-fin, one specifies the temperature, T = 


Ts». It is then necessary to estimate B/2, aided by the ap- 
proximate solution of Section 4 (e.g., assume a value for B/2 about 
10 per cent greater than that given by the approximate solution 
for the same 7:,,). If as a result of the computations the tem- 
perature at zr = D/2 is Tim, then the assumption of B/2 is con- 
firmed; if not, then the computations must be repeated using an 


dT’ 
adjusted estimate of B/2. The magnitude of (*) will 
dx z=D/2 


result from the computations, and the heat rejected from the fin 


(=) 
dz r=D/2 


We will examine the procedure by considering the example of 
section 4 (a fluid condensing at 1800 deg R). 7, = 1800 R, let 
T,; — T, = 200 R; therefore T, = 1600 R. Assume B/2 is 10 
per cent greater than the approximate solution yields (Table 3), 


per unit length, at 7,,, is (Q/L), = —kt 


ie., B/2 = 1.1 & 0.405 = 0.45 in. The fin width B/2 is broken 
into five intervals, Ar = 0.90in. Starting at mid-fin, — apt d 


= 0, T = 1600 and F,, can be obtained by replotting the informa- 
9 
2. — = 18; ie., a curve of 
0.5 

F, versus z along the fin. 7” is computed from equation (13), 
and —AT7’ = 7” Ar except for the first and last steps, for which 
—AT’ = '/,T" Az, to establish the slopes at the mid-poirt of each 
increment; —7’ = 2(—AT’). AT = —T’Az and the tempera- 
ture at the end of the first increment is the temperature at the 
beginning plus A7’. The summary for this example appears in 
Table 4, where 7, is computed to be 1804.65 deg R, which com- 
pares closely with the original value of 1800 deg R. (Q/L),(exact) 
is compared with the approximate value for a range of B/D in Fig. 
10. The product F,7* is shown in Fig. 11; this represents the 
effectiveness of fin at any point. Note the adverse effect of poor 
view factor near the tube-fin junction. Fig. 12 compares the fin 
temperature distribution for the approximate and exact solution. 
Fig. 13 compares the maximum heat rejection per unit mass for 
the exact and the approximate solution. Note that the optimum 
B/D ratio shifts slightly. 


tion given in Fig. 17 for B/D = 
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W. Optimization of a Conventional-Fuel- 
Fired Thermoelectric Generator 


E. V. SOMERS? 


Introduction 


THERMOELECTRIC GENERATOR involves the direct 
production of electrical power by flow of heat into and through a 
thermocouple system. Thus, as depicted in Fig. 1, if a series 
circuit of thermocouples of alternate materials ‘“X”’ and “Y”’ 
are placed between a heat source and a heat sink, heat will flow 
into the hot junctions and part of it will be converted to electric 
power, as indicated by the + and — signs at the electric ter- 
minals, and part of it will be rejected to the sink at the cold june- 
tions. 

Shown in Fig. 2 is a schematic diagram of a cross section of a 
thermoelectric generator duct. The diagram illustrates the ex- 
tension of the simple scheme of Fig. 1, with the heat source that 
of the heated gas flow rate in the center of the duct and the heat 
sink that of the coolant flow rate in the coolant jacket. The 
coolant is assumed to remain at ambient temperature along the 
length of the duet. However, the gas drops in temperature as it 
flows through the duct because of the heat rate transferred by 
convection to the thermoelectric elements. The convected heat 
rate at any point along the length of the duct can be computed 


'Heat Transfer Section, Mechanics Department, Westinghouse 
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? Manager, Heat Transfer Section, Mechanics Department, West- 
inghouse Research Laboratories. 

Contributed by the Heat Transfer Division and presented at the 
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Nore: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those of 
the Society. Manuscript received at ASME Headquarters, March 2, 
1959. Paper No. 59—SA-14. 


Nomenclature 


in terms of the over-all transmission coefficient and the gas-to- 
ambient temperature drop. The over-all transmission coefficient 
is a function of the gas and coolant convection coefficients, ap- 
proximately constant along the length of the duct, and the local 
thermal conductance of the thermoelectric elements. Since 
the local thermal conductance is a design variable dependenton the 
local element thickness, the design of the element thickness along 
the length of the duct will determine (1) the local gas temperature 
and the local convected heat rate; and (2) the local thermoelectric 
power produced. Thus, for a given duct geometry and given gas 
conditions of inlet temperature and mass flow rate, it is desirable 
to determine what temperature distribution of the gas along the 
length of the duct will optimize the total thermoelectric power 
produced and what distribution of thermoelectric element thick- 
ness is necessary to obtain this gas-temperature distribution. 


i Heat Source | 
i 


Flow 


*Hot” 
Junctions 


“Cold” 
Junctions 


cross-sectional area normal to 


heat flux rate in thermoelec- 
tric elements, ft? 

equations (15) 

specific heat of heated gas at 
constant pressure, Btu/lb 
deg F 

duct diameter, ft 

arbitrary function of gas tem- 
perature and its derivative 
with respect to duct length, 
equation (12) 

enthalpy rate of heated gas, 
Btu/hr 

gas heat-transfer coefficient, 
Btu/hr ft? deg F 

heat-transfer coefficient to sink 
coolant at an ambient tem- 
perature, Btu/hr ft? deg F 

kinetic energy rate of heated 
gas, Btu/hr 

thermoelectric element thick- 
ness, ft 

thermoelectric element length 
optimizing power, ft 
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duct length, ft 

mass flow rate, lb/hr 

potential energy rate of heated 
gas, Btu/hr 

thermoelectric power, Btu/hr 

optimum thermoelectric power, 
Btu/hr 

heat flow rate, Btu/hr 

heat rate added from fuel to gas 
flow at duct entrance, Btu 
hr 

temperature of heated gas, deg 
R 

hot-junction temperature of 
thermoelectric elements, deg 
R 

cold-junction temperature of 
thermoelectric elements, deg 
R 

ambient temperature of sur- 
roundings, deg R 

inlet bulk gas temperature, deg 
R 


exit bulk gas temperature, deg 
R 


over-all transmission coefficient 
Btu/hr ft? deg F 


distance co-ordinate along 
length of duct, ft 


linear gas temperature gradient 
in the duct, deg F/ft 

linear gas temperature gradient 
optimizing power, deg F /ft 

Carnot efficiency 

thermocouple efficiency 

modified thermocouple effi- 
ciency, 7 , deg R 

thermal conductivity of ther- 
moelectric material, Btu/hr 
ft deg F 

thermal efficiency 

optimum thermal efficiency 
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Fig.2 Cross secti 


of thermoelectric generator duct 


Determination of Electric Power in Terms of Design 
Parameters 


The design parameters of the thermoelectric duct are those of 
geometry, duct diameter D, and duct length L; material, thermo- 
electric parameter €.¢.*, and specific heat of gas c,; and operating 
conditions, inlet gas temperature 7, gas flow rate M, coolant 
temperature 7',, and convection coefficients h, (gas) and h, 
The local thickness of the thermoelectric elements 
along the length of the duct determines the local differential heat 
rate dQ and the loeal hot and cold junction temperatures 7’, and 
7, and these in turn determine the local differential thermoelec- 
tric power dP. The integration of dP along the length of the 
duct yields the total eleetrie power P. 


(coolant ). 


Let P denote thermoelectric power; in differential form [1 


dP = dQ (1) 


where ¢€, is the Carnot efficiency (7, — 7T.)/T, and €<¢. is the 
thermccouple efficiency. 

From the first law of thermodynamics, the energy content of 
the heated gas flowing in the thermoelectric duct is described by 


dQ — dW = dH +dKE + dPE 
dW =dKE = dPE~0O 


dQ = dH = (2 
or 
dQ = Mc, 
“s Or 


The heat-flow rate described by equation (2) equa's that heat 
rate transferred radially through the thermoelectric elements 

The differential thermoelectric power, dP, may be rewritten 
from equation (1) as, 


Pe aQ (3) 
T, 
By definition‘ let 
€vc. 
dP = — (4) 


* Numbers in brackets designate References at end of paper. 

‘ This definition is approximately correct (very close) for thermal 
efficiencies around 10 per cent; it results from a binomial expansion 
of loffe’s value of M [2] with loffe’s value of z/4 equaling the value of 
*. is a material parameter. 


€t.e. 
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For any section of thermoelectric elements, the Fourier law for 
heat conduction states that, as an approximation (for small 
values of the €,€:... product), [3] 


T, —f. 
dQ = dA (5) 


Therefore, an alternate expression exists for dP 
dP = &...* : dQ? (6) 
dA 
A further representation of dQ is given by 


(T, — 


dQ = U(T, — T,A = (7) 


from which 


dP = \ UF ue dA (8) 


An auxiliary condition given by equations (2) and (7), withdA = 
is 


Me dT, — T.) (9) 


which can be solved for (2) 


— T,) 1 
Ur) = + ) 
Me 


(10) 
dr 


This expression for [(z) can be substituted into equation (8) for 
dP. 


L 
iT 
P= )-Me,—(T, T 
1 | 
dT ,\? 7 


dx 


Henceforth, € 
of the duct.® 


* will be assumed to be constant along the length 


Optimization of Electric Power with the Optimum Tempera- 
ture Distribution of the Gas 


The problem now is to maximize the integrated expression for 


P = F\T,, with T, at (12) 
0 dx 


and 7, unspecified at c = L. This requires that F (7, ‘) 


satisfy the Euler-Lagrange differential equation [1] 


5 The value of &..* can be assumed constant using temperature- 
averaged values for the thermal-conductivity-electrical-resistivity 
product and the mean value of the Seebeck coefficient at the hot and 
cold junction temperatures [4]. Although the hot and cold junction 
temperatures vary in this problem, an a priori estimate of &...* equal 
to a constant introduces little error. 
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oF d ( or ) 0 
oT, dz \oT,,’ 
and, in this case, the end point condition that [5| 


——-=Q@ at 


r=L. 


1 \ (Mc,)* 
~) 7D 


L 
P= {A,\T,(T, — Te) + dx | 
| 
| 


oF 
oT,’ 


= A(T, + 2A;T,’ 


By substituting into equation (13), it can be shown that 
T, = To + Br (17) 


The end point condition, equation (14), now determines the op- 
timum values of 8, 8*, from equations (15) and (16) 


(To 
2Mc 1 
L P 
(; 


= 7) + Br, equation (11) becomes 


B* = — (18) 


with 7, 
L 

0 


dr (19) 
This integrates to 


L? 
*Mc,B( To — T,)L — €.c.*Mc,8* = 
(A 3 } 
h, 


*Me, [To — — T, 


Mc, ( 1 


—(T,) — + 


Substitution of 8 = 8*, yields 


P* = *Me,( = )? 


, L) 
h, 


The expression for P* can also be written in terms of 7',*. 
ease L is specified and 7,,* is obtained from the equation 


In this 
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2Mc, fh 
L+—{— 


2A 
The solution for L + a ( 2. + a substituted into equa- 
mD \h, h, 


tion (21) gives 


— T.XTo — T.*) 
2 


= Etc. (23) 


It is to be noted that satisfaction of the Euler-Lagrange dif- 
ferential equation is a necessary but not sufficient condition for 
maximizing the power integral. By using exponential and cosine 
gas temperature distributions for comparison, it can be shown 
that 7, = 7) + Bx does not minimize the power. 


Determination of Optimum Thermoelectric Element 
Thickness 


In the optimization of the thermoelectric power, the condition 
used was that 


rD(T, — 


24 
( 1 l 1 (24) 
h, h, 


This equation can be solved for I(r) 


ArD(T, — 
dT, 


Uz) = 


For T, = To + Bx 
_ AWD((To — T.) + Br) _ 


Ur) = 


Substitution of 8* for 8, per equation (22), yields 


Iz) « Ar@D(L — xr) 
Me, 


Thermal Efficiency Under Optimum Conditions 
The thermal efficiency of interest is defined as the quotient of 
the optimum thermoelectric power by the heat-addition rate to 
the system. Thus, for fixed values of M, é..*, To — 7, and L, 
Qe 


By equation (23) 


— 


9 


— T,*) 


pP* = 


and since Qi, = Me,(T> — T,) 
€* = To — T,) 

As the optimum duct length approaches infinity, the exit tem- 
perature approaches ambient. In this case the efficiency ap- 
proaches its maximum value 
(31) 


e*(max) = '/2€t¢.*(T> — 


For the case in which the exit temperature approaches the inlet 
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4 
oF 
(14) 
4 * ( 1 (15) 
h, 
Then, 
dT 
Mc, —! om 
dx 
: 
P 
dx. 
j 
(2 
+ 26) 
h, h, 
: 
wt 
} 
+ | _ 99 
(T, — T, | 
(21) 


temperature the optimum duct length and the efficiency €* both 
approach zero. 


Summary 


The purpose of this report has been to investigate conditions 
for optimizing the thermal efficiency of the design of a conven- 
tional-fuel-fired thermoelectric generator. Equations have been 
developed for the following: 


1 The temperature distribution of the heated gas in a circular 
duct type of generator, that maximizes the thermoelectric 
power output, with respect to the temperature distribution within 
the duct and the exit gas temperature at the end of the duct. 
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2 The thermoelectric power developed, based upon these 
optimizing temperature criteria. 

3 The required variation of thermoelectric element thickness 
required to satisfy both optimizing temperature criteria. 

4 Thermal efficiency under optimum conditions. 
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Heat Transfer From a Rotating Disk 
to Fluids of Any Prandtl Number 


E. M. SPARROW? and J. L. GREGG’ 


Introduction 

THE LAMINAR FLOW about a rotating disk situated in a 
large body of quiescent fluid was first analyzed by von Karman 
{1].2. The heat transfer from such a rotating disk has been 
studied by Millsaps and Pohlhausen [2] for fluids which have 
Prandtl numbers in the range 0.5 < (c,/c,)Pr < 10. For gases 
(Pr = 0.72), the effects of compressibility were examined by 
Ostrach and Thornton [3]. In the present investigation, the re- 
striction on Prandtl number is lifted, and heat-transfer results are 
obtained for fluids of all Prandtl numbers. 


Analysis 


The heat-transfer process in the fluid is governed by the con- 
servation of energy principle, which takes the following form in 
cylindrical co-ordinates 


oT Vs OT o7 
— kV?7T 1 
( + + V, ~~) (1) 


where T is the static temperature and Y? is the Laplace operator. 
The velocities may be rephrased in terms of von Karman’s simi- 


larity variables as follows: 
1/2 
rwF(n), Ve = rwG(n), 


= (wy 


‘2H(n) (2a) 
and in addition, a dimensionless temperature may be defined as 
-—T. (2b) 


where 7, and T,, respectively, represent the surface temperature 
(a constant) and the ambient temperature. Introducing these new 
variables into the energy equation, there is obtained 


6” = (Pr)H6’ (3) 


The primes denote differentiation with respect to 7 and Pr is the 
Prandtl number of the fluid. From the definition of @ as given by 
equation (2b), it is clear that the boundary conditions are 


a0) =1, 0+0 as n> & (4) 


Numerical solutions of equation (3) subject to the boundary 
conditions (4) have been carried out on an IBM 653 electronic 
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computer for Prandtl numbers of 0.01, 0.1, 1, 10, and 100.4 The 
temperature distributions corresponding to these solutions are 
presented in Figs. 1 and 2, respectively, for the low and high 
Prandtl number ranges. Also shown on the graphs are the values 
of the velocity function H used as input data in the solution of 
equation (3). Inspection of these figures reveals that, for low 
Prandtl] numbers, the thermal boundary layer is much thicker 
than the velocity boundary layer. The opposite characteristic is 
displayed for high Prandtl numbers. This suggests a procedure 
for obtaining asymptotic solutions. 

First, as the Prandt] number becomes very small, Fig. 1 indi- 
cates that H is essentially constant throughout the thermal 
boundary layer. Under these circumstances, the solution of equa- 
tion (3) is 


= ePr Hl(@)n (5) 


where = 


—0.88447. In particular, the temperature 
derivative at the 


wall, needed in computing the heat transfer, is 


(d0/dn)q-0 = —0.88447 Pr (5a) 


From Fig. 2, it is seen that, as the Prandtl number approaches 
very large v toe the thermal boundary layer is confined to a 
smaller and smaller portion of the velocity boundary layer. This 
prompts us to write H in terms of a series expansion about 9 = 0. 


So, 


H"(0) 
9 


H = H(O) + H'(O)n + —_?+... 


Since Ho) = H'(0) = 0, we Ww rite 


2 
H = H"(0) 7 


for small values of 7. A solution of equation (3) corresponding 
to this approximation for H is 


[—Pr H”:0) 
f, nexp 
(6) 
f, exp (— 


The dimensionless temperature derivative at the wall associated 
with this solution is 


O=1- 


= H"(0)/6)'*/T(4/3) = 0.62048 Pr'/* (6a) 


where the value — 1.02046 has been used for H"(0), and repre- 
sents the gamma funetion as found in numerous mathematical 
tables. 


Heat-Transfer Results 
The local rate of heat transfer from the disk to the fluid may be 


evaluated using Fourier’s law, i.e., 


q = —k(OT/dz) (7) 


*The values of H(n) needed to carry out the solutions were ob- 
tained by re-solving von Karman's differential equations. 
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Fig. 1 Temperature distributions for low Prandtl number fluids |H/H( ~) 
is input function for temperature solutions) 


Introducing a local heat-transfer coefficient h by the definition 


h = q/(T, — T.) 


and evaluating equation (7) in terms of the variables of the analy- 
sis, there is obtained 


(7a) 


( dé 
dn nm 
Since (d0/dy)y.» depends only on the Prandtl number, it is seen 


that, for a given fluid, the variation of the heat-transfer coefficient 
with angular velocity is given by 


h~o'”? 


Further, h is constant over the disk surface. 
Using the numerical solutions of equation (3), the dimension- 
less heat-transfer results are listed in Table 1. 


Table 1 Dimensionless heat-transfer results 


0.01 00087051 

0.1 0.076581 

0.39625 
10 1.1341 
100 2.6871 


For the limiting situations of very low and of very high Prandtl 
numbers, the following asymptotic expressions for the heat trans- 
fer may be written utilizing equations (5a) and (6a) 


= 0.88447 Pr, (8a) 
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Fig. 2 Temperature distributions for high Prandtl number filvids |H/H( 
is input function for temperature solutions 


| 


onless heot transfer 


Dimens 


Heat-transfer results 


0.62048 Pr— (8b) 

A graphical presentation of the heat-transfer results is given in 
Fig. 3. The results have been plotted in two different ways, one 
appropriate to low Prandtl number fluids and the other appro- 
priate to high Prandtl number fluids. 
also shown (dashed). 


The asymptotic lines are 
It is seen that the low Prandtl number 
asymptote already closely coincides with the computed curve at 
Pr = The high 
Prandtl number heat-transfer results approach their asymptote 
somewhat less rapidly, the difference between the computed and 
asymptotic curves being about 7 per cent at Pr = 100; this 
deviation is not large from the practical point of view. The sig- 
nificant fact is that the asymptotic curves provide sufficiently 
accurate heat-transfer results outside the Prandtl number range 
for which the numerical solutions have been obtained. 

By utilizing either the numerically computed results or the 
asymptotic expressions, heat-transfer calculations can be carried 
out for fluids of any Prandtl number. 


0.01, the deviation being only 1.5 per cent. 
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An Approximate Treatment of Unsteady 
Heat Conduction in Semi-Infinite Solids 
With Variable Thermal Properties 


KWANG-TZU YANG' and ALBIN SZEWCZYK? 


THE PROBLEM of property variations in unsteady one-di- 
mensional conduction analysis has become increasingly impor- 
tant in practice, especially in view of extreme temperature 
changes encountered in many applications. As the boundary 
temperature conditions are usually arbitrary, it is therefore de 
sirable to have a relatively simple, yet accurate, calculation pro- 
cedure for such a general problem so that the variation of proper- 
ties may readily be taken into account. The primary purpose of 
this note is to present such a calculation procedure. The in- 
tegral method is well known in boundary-layer theory; its ex- 
tension to conduction analysis, however, has only attracted at- 
tention recently. Goodman [{1],* utilizing temperature profiles 
approximated by polynomials, solved integral equations involv- 
ing a phase change. Very recently, Reynolds and Dolton [2] 
presented a more general analysis applicable to a variety of con- 
duction problems, also using temperature polynomials. The 
present procedure may be considered as a significant improve- 
ment over that of Reynolds and Dolton [2] for treating conduc- 
tion in semi-infinite solids with arbitrarily varying surface tem- 
peratures, since their method fails to yield physically sensible 
results for large dependence of thermal properties on tempera- 
ture. 


By introducing the following dimensionless quantities: 2* = 
t* = at/L?, T* = T/T ;, k* = k/ko, C* = C/Cu, a = ko/Co; 
and 6 = (7 — T,;)/(T, — T;), where x is the depth into the solid, 
T the time variable, L a characteristic length, & the conductivity, 
C the volumetric specific heat, and subscripts 7, s, 0, indicate 
initial, surface, and reference conditions, respectively, the usual 
one-dimensional unsteady conduction equation with variable 
properties may now be written as 


= (T,* — 1) ‘G -) (1) 
or* 


which is to be solved with the conditions 6(r*, 0) = Ofor0 < a* 
©, 0, t*) = 1 for t* > 0, and & ~, = Ofor’(* > 0. Now 

' Department of Mechanical Engineering, Notre 
Dame, Notre Dame, Ind. Assoc. Mem. ASME. 

‘ Now at Institute for Fluid Dynamics and Applied Mathematies, 
University of Maryland, College Park, Md. Assoc. Mem. ASME. 

* Numbers in brackets designate References at end of paper. 

Contributed by the Heat Transfer Division of THe AMERICAN 
OF MECHANICAL ENGINEERS. 

Nore: Statements and opinions advanced in papers are to be 
understood as individual expressions of their authors and not those 
of the Society. Manuscript received at ASME Headquarters, March 
2, 1959. 


University of 


Journal of Heat Transfer 


integrating Equation (1) with respect to x* from 0 to © and im- 
posing the condition C* = 1 + B7T* = 1 + B[(T,* — 1)6 + 1], 
where @ is a constant, we obtain: 


) + + B(T,* — 1)6.*] 


B dé.* 
1 = = — 


@ 


where 6,* = f, 6dr* and 6.* = f, 6%dr*. In order to solve 


Equation (2), it is necessary to introduce an assumed unsteady 
temperature profile. Because of the limitation on using poly- 
nominals as indicated earlier, the present analysis utilizes an ex- 
ponential form of the type 


yee z (3) 


where 6*((*) and J(¢t*) are two unknown functions. Physical 
consideration requires J to assume negative values, and in fact 
it ranges only from —1 to 0, the former being corresponding to 
no heat flux at surface. It is also seen that usual asymptotic 
conditions as x* — © are automatically satisfied in the assumed 
profile. With Equation (3) substituted in Equation (1) for z* = 0 
and in Equation (2), and 6*(¢*) eliminated, there results the final 
integral equation: 


Ox*, t*) = 


dJ P(J, t*) — M(J, t*) 
dt* NJ, t*) 


where 
P(J, t*) = k,*Fs 
M(J, t*) - ( + G*F.) 


1 
+ 4 (E*'F, + G*'F.){B[2 + F(T,* — 1)] + 2} 


N(J, t*) = (E*F,’ + G*F,')} B(2 + FA T,* — 1)] + 2} 


+ B — 1)F (ETF, + s) 


G*(t*) = 


(1 — J?)? 


1) 


dt* 
= (2 + 1) J)? FAJ) 
J? +2 
F(J) = 


FJ) =1—- J? 


and the prime indicates derivatives. A limiting process of Equa- 
tion (4) as ¢* — 0 yields the initial condition /J(0), which is de- 


termined from 


1 E*’ 
F;[J(0)] = + + B ( ) (5) 
2 


Equation (4) may now be solved by any convenient numerical 
procedure, yielding J(t*). Finally, unsteady temperature pro- 
files are given by Equation (3) with 6*(*) = (E*F, + G*F,)'*/(1 
— J) and surface heat flux q, by (q,L/keT,) = (T,* — 1k,*C + 
J)/6*. 
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Fig. 1 T ient temperature distribution with stepwise surface tem- 
perature change and variable thermal conductivity 


— variable properties 
-—— constant 


30 


perature distribution for the numerical example 


For a demonstration of the accuracy, the present procedure is 
compared with an exact solution given by Yang [3] for the 
problem of a semi-infinite solid undergoing a stepwise surface 
temperature change with 8 = O and k* = 1 + a@@. Such a com- 
parison is shown in Fig. | fora = +0.5, 0, and —0.5, where y = 
factory, especially in the region where large temperature varia- 
tions occur. Incidentally, during the course of this study, a 


It is seen that the present procedure is entirely satis- 
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variable properties 
constant 


Fig. 3 Surface heat-flux characteristics for the numerical example 


fourth-degree polynomial was also tried, and it was found that no 
physical’y sensible results were obtainable for this problem with 
a = +0.5. Finally, to show an application of the present proce- 
dure to a general problem, an example with 7,* = e”, 8 = 0.034, 
and k* = 1 — 0.0577T* is chosen. The results are shown in 
Figs. 2 and 3. Also shown are the corresponding curves for con- 
stant properties. It is interesting to note that in this particular 
problem, the unsteady temperature distribution is affected by the 
property variation more than the surface heat flux is. It is be- 
lieved that the present procedure may also be extended to more 
general conduction problems. 
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